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ABSTRACT 


This thesis involves investigation of linear filtering 
models as a means of generating texture in images. Various 
autoregressive filter models are used to generate various 
textures, and the results are analyzed to determine rela- 
tionships between filter parameters and texture characteris- 
acs . A two-dimensional counterpart to the autoregressive 
integrated moving average (ARIMA) model from one-dimensional 
time series analysis theory is developed and tested for 
texture modeling applications. All these models are driven 
by white noise, and to the extent that real images can be 
reproduced this way, advantages in image texture transmis-— 
Sion could be realized. Results of this work indicate that 
the purely autoregressive models work well for some types of 
image textures, but that for the textures studied the ARIMA 


model is not particularly suitable. 
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The purpose of this thesis is to investigate the types 
and quantity of image textures generated using a_ two- 
dimensional (2-D) extension of the Autoregressive Integrated 
Moving Average (ARIMA) model. For the one-dimensional 
(e.g., time series) case, the theories and formulas describ- 
ing this model are outlined in Box and Jenkins [Ref. 1:pp. 
85-103]. The 1-D ARIMA model is useful when the time series 
to be modeled is not stationary but exhibits some 
homogeneity in the sense that, except for statistical 
differences between parts of the time series, these 
different parts of the process behave similarly. In these 
cases some suitable difference of the process may be 
stationary, and hence may be accurately modeled by an 
Autoregressive Moving Average (ARMA) OF a purely 
Autoregressive (AR) model. The resulting stationary time 
series (generated by an appropriate ARMA or AR filter with 
white noise input) is applied to an integra-tion or 
Summation filter (the inverse of the difference operation) 
to generate the original nonstationary time series. (Ref. 
1l:p. 85] Figure 1-1 shows a block diagram of this process. 

This work attempts to extend these concepts to two- 
dimensional signal processing. In order to simplify the 


model, the moving average (MA) portion of it will be 
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Figure 1-1 Block Diagram for the Autoregressive 
Integrated Moving Average Model 


eliminated (i.e., no zeros in the filter Z transform) so 
that only purely AR models will be considered for stationary 
image generation. The procedures for modeling image 
textures using AR models with white noise input are well 
established [Ref. 2:pp. 454-456]. However, a suitable two- 
dimensional difference operation and its inverse must be 
found to implement the concepts outlined above. 
The research is divided into four areas: 


1) Investigation of the various types of image textures 
generated using AR models where filter coefficients 
and size are determined a) arbitrarily, b) using a 
two-pole separable model, and c) using a four-pole 
separable model. Separability refers to the fact that 
the Z transform of the AR filter can be factored into 
components representing each dimension or direction of 
the image. 


2) Selection of a difference operator and a realizable 
inverse (integration or summation) filter. 


3) Application of the above autoregressively generated 
images to the summation filter, and evaluation of 
these results. 


4) Attempted reproduction of actual images textures using 
AR models whose coefficients are determined using the 
statistics of the image, and comparison of these 
results to those obtained by a) applying the 
difference operator to the real image, b) finding the 
coefficients of the AR model that reproduce the 
difference image, and c) applying the difference image 
to the summation filter. This comparison was intended 


to discover what improvements, if any, may be realized 
using the ARIMA model vice a purely ARMA (or AR) 
model. 

The remainder of the thesis 1S organized as follows. 
Chapter II contains methods and results of investigating 
various autoregressive image models. Chapter III deals with 
the formulation, development, and testing of the two- 
dimensional summation filter. Chapter IV contains the 
results of applying various AR-generated images to the 
summation filter and addresses the application of the ARIMA 
model to real image textures. Chapter V outlines conclu- 
sions on the results and applicability of the ARIMA model. 
Although the ARIMA modeling was not highly successful in 
reproducing the textures studied here, plausible reasons are 
given for their failure and conjectures are made about those 
circumstances where the model would be more successful. 
Appendix A provides information on the computer algorithms 
used to implement the equations governing the above 
processes. Appendices B through G contain derivations of 
spectral and autocorrelation equations, and the correspond- 
ing graphical results, governing the AR processes in Chapter 
Tl. Appendices H and I contain graphical results associated 
with the inverse filter development in Chapter III. 

Image data were generated using computer programs 
written in FORTRAN, compiled using Version 4.5 under the 
VAX/VMS Version 4.4 operating system. The images were 


displayed on the COMTAL Vision One/20. The gray level 


intensity range of pixel values is 0O (darkest) to 255 
(lightest), so the image data generated had to be scaled to 


that range for display (see Appendix A). 


0 


Il. THE AUTOREGRESSIVE IMAGE MODEL 


A two-dimensional signal (such as an image texture) can 
be modeled using a two-dimensional AR model with white noise 
input. The governing equations in the spatial domain are of 
mie following form ~Ref. 3:pp. 325-326]: 

P= O=v 

y(n,m) = - ) ) aijy(n-i,m-j) + w(n,m) (2-1) 
i=0 43=0 
(1,j)4(0,0) 
where y(n,m) iS a Signal representing the generated image 
texture at pixel location (n,m), aig is the filter coeffi- 
cient matrix, and w(n,m) is a two-dimensional white noise 
Signal. The system function corresponding to the filter of 
Boa. (2.1) is given by 
iy. = : 1 
pene 72) = Dh 92 Hp 2p gg 2] ote pag 4 Za ) 
: W(Z, 12.) 





(2) 


where Z; and Zz are the Z transform variables corresponding 
to spatial coordinates n and Mm. Ideal white noise has an 
autocorrelation function that is an impulse and a flat (con- 
stant) power spectrum with magnitude corresponding to the 
variance of the white noise process [Ref. 4:pp. 22-26]. 
Therefore determination of the filter coefficients aj4 will 


define the generated image process. Procedures will be 


ia 


outlined later to estimate the coefficients from real image 
data. At this point analytical methods will be used to 
select these coefficients and the resulting images will be 
studied. Figure 2-1 shows an example of a white noise input 


image. 
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Figure 2-1 White Noise Image 


A. ARBITRARILY SELECTED FILTER COEFFICIENTS; P = 2) 508—a2 
In order to get an initial idea of what types of images 
might be generated using a 2 x2 AR filter with white noise 
input, filter coefficients were at first selected arbitrari- 
ly, but subject to a stability constraint. The primary con- 
straint on coefficient selection is that of filter 
stability. Using the DeCarlo-Strintzis Theorem dealing with 
multidimensional filter stability [Ref. 3:pp. 197-198], 
alternately setting z,; = 1 and Zz = 1 and determining the 
location of the pole in the remaining dimension will 


indicate whether or not the filter is stable. If the 


ae2 





magnitude of the pole in the remaining dimension is less 
than 1, the filter is stable. Even with this condition, 
however, there are an infinite number of possible filter 
coefficient combinations. The additional constraint of 
€a10 201 can be used, and comparisons of results using 
various values of aj, can be made. 

Figure 2-2 shows the form and directionality convention 
used for the autoregressive filter, along with the corres- 


ponding difference equation and its Z transform. 


h(n,m) 
m(Z2) aq] 211 I L 
l - ; y(n,m) = - ) L Se cue + w(n,m) 
an) —() j= 
nee) ) (1,3)4#(0,0) 
iu 


¥ (24,25) = H(Z),Z5)-W(Z),Z5) = jay OC OWl2y/2>) 


-l -l 
l+a,.Z. +a,,.Z. +a 1 Z5 


1071) 01-2 ae 


Figure 2-2 Autoregressive Filter Impulse Response, 
Difference Equation, and Z Transform 


Although it is difficult to make precise predictions in 
two dimensions, one can expect that the sign and magnitude 


Of a}9 OF ap, would influence the correlation between pixels 


3 


in the corresponding directions. For example, a positive 
value for aj, 9 might be expected to yield an image with 
substantial variation in the n direction (low correlation), 
particularly if the magnitude of ajo is near 1. A negative 
value for ap, with magnitude near 1 might yield an image 
with lower variations in pixel intensity (high correlation) 
in the m “darectione Since the filter is not necessarily 
separable (i.e., the denominator of H(z2,,2Z2) cannot be 
factored into the form D(Z1)°D(Z9)), conclusions drawn from 
this line of reasoning may not be completely correct. 
Initial attempts at generating images with arbitrarily 
selected coefficients yielded rather uninteresting results 
having very little contrast or discernible pattern. 
Continued experimentation with combinations where 
€a190 = ao, < O eventually yielded more interesting image tex- 
tures. Figures 2-3 and 2-4 show the results of using the 
constraint a19 = 491 = -0.35 and various values of aj, for 
the filter coefficients. For positive values of aj, the 


images are rather "grainy," with higher magnitudes yielding 


a somewhat "finer" graininess. There are also some overall 
intensity differences observed. For the negative values of 
€ai1, the results are much more interesting. As the magni- 


tude increases, there is a gradually more noticeable upper 
left to lower right orientation of the image texture, and 
the variations from lower left to upper right become 


smoother as well. Using aj9 =& 1 = -0.38 and a4 4 = -0.24 


14 





Figure 2-3 Images Generated Using Arbitrarily 
Selected Filter Coefficients 
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aq1 = -0.3 


Figure 2-4 Images Generated Using Arbitrarily 
Selected Filter Coefficients 
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yields minor variations in texture pattern and overall image 
intensity when compared to the previous case. 

Using 419 = a9, = 0.35 and varying aj, from 0.0 to 0.8 
(Figure 2-5), the images obtained deviate very little from 
the mean intensity value, and possess minor differences in 
graininess. With these positive coefficients some negative 
correlation might be expected, and the fact that these 
images are "grainy" indicates the existence of some negative 
correlation or high spatial frequency characteristics. On 
initial examination, however, the low contrast of the 


generated images tends to obscure the observed graininess. 


B. TWO POLE, SEPARABLE AUTOREGRESSIVE MODEL 

In general, it is difficult to relate the nature or 
properties of a two-dimensional filter to the precise nature 
of an image texture that may be generated when white noise 
is applied to that filter. In order to simplify the effort 
and to obtain a better understanding of the problem, the 


case where the filter (and resulting image texture) are 


separable is considered. For the two pole separable case 
considered here, the filter transfer function can be 
factored into expressions in 2, alone and Zz alone. The 


Expressions in 2, and Z5 each have one pole on the real axis 
in their respective Z domains. Figure 2-6 illustrates the 
filter structure, the corresponding difference equation, and 


mes 4 transform. 
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Figure 2-6 Autoregressive Filter Form, Difference 
Equation, and Z Transform 


Here it is relatively easy to relate stability of the 
filter to the location of the poles in the zy, and z»y planes. 
Since the quarter plane filter is separable and the 
components are causal, one-dimensional filter stability 
theory can be used to state that the poles in each plane 
must have magnitude less than 1 to ensure filter stability. 
Figures 2-7 through 2-9 show images resulting from this 


model for various values of a}39 and ap,. Note that the sign 


convention used for the filter difference equation and Z 
transform results in poles on the negative side of the real 
axis for positive values of a)9_0r a9;. and vice versa. 
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Figure 2-7 Images Generated Using a Two-Pole 
Autoregressive Model 
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Figure 2-8 Images Generated Using a Two-Pole 
Autoregressive Model 
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Figure 2-9 Images Generated Using a Two-Pole 
Autoregressive Model 


Be 


Careful comparison of the images resulting from various com- 


binations of aj19 and ao, =\XIleads to the _ following 


observations: 


1) 


2) 


3) 


4) 


3) 


6) 


When poles are located in the same place in the 2, and 
ee eee se ie ec cen ee Sie OL st hie wu G@ealmecax 1S), 
magnitudes near 1 yield a fine graininess with patchy 
areas and low overall contrast. As the magnitude of 
the pole decreases, the graininess becomes more coarse 
and the result is more like the original white noise 
input. No directional quality in the image pattern is 
observed. 


When poles are located in the same place on the posi- 
tive side of the z, and Zz real axes, a somewhat dif- 
ferent result is observed. For magnitudes near 1, an 
image of patchy light and dark areas results, with 
differing amounts of correlation between pixels in 
different areas. Slightly discernible "lines" in both 
the horizontal and vertical directions are also 
observed. For lower pole magnitudes on the positive 
Side of the real axis, the decreased effect of the 
filter on the white noise input is again observed. 
This result is more like the white noise and has more 
contrast than the corresponding result using poles on 
the negative side of the real axis. 


For pole placements in the z,1 and zZ> planes which are 
on the negative side of the real axis and are of 
unequal magnitude, the results have a very fine 
graininess and low contrast. Some slight direction- 
ality is observable in the image patterns, with lower 
frequency variations evident in the direction corres- 
ponding to the pole with smaller magnitude. 


For pole placements in the Z, and Z> planes which are 
on the positive side of the real axis and are of 
unequal magnitude, much more directionality and varia- 
tion is observable in the image pattern. 


As the poles are placed on opposite sides of the real 
axis and are separated by a greater distance, direc- 
tionality becomes more evident (with higher frequency 
variations in the direction of the more negative 
pole). As the pole separation becomes greater and as 
the pole magnitudes become closer to 1, smoother sinu- 
soidal variations are evident. 


When the pole magnitudes are equal and have opposite 
sign, the image generated using pole magnitudes close 


De 


to 1 exhibits high frequency sinusoidal variations in 

the direction of the negative pole. The image gener- 

ated with the lower magnitude poles, as would be 

expected from the above results, resembled the 
unfiltered white noise. 

In order to explain these image patterns analytically, 

analysis of the power spectrum and autocorrelation function 

of this process is useful. Since this model is separable, 


the analysis can be conducted in each direction separately. 


The power spectrum is defined by [Ref. 4:pp. 24-34): 
Sy(w) = o2|H(eH)|2 = ¢ 2H(eI4H(e7I%) (2. 3a) 
where 
H(eJ¥) = H(z) [za (2.3b) 


and for this case 


H(z) = — (2. Sip 


1 +a0z 
Here o% is the magnitude of the white noise power spectrum. 
We can assume that o% = 1 with no loss of generality of the 
results, since o? does not affect the Shape of the frequency 
response. 
The autocorrelation function is related to the filter 
transfer function through the equations [Ref. S:pp. 391- 


SOE: 


24 


ney ee ee il (02) 0) (2.4a) 
SC emo Osy artey, (2) 
Ry(2) = o° a Euan) eta (ines) (2.4c) 


Specific forms of the power spectrum and autocorrelation 
function are given in Appendix B. Since Ry(&) = Ry (-* ) 
(Ref. 5:p. 388], calculating the expression for Ry (2), n=O 


is not necessary. From Appendix B, the results are: 


Sy (W) = : D) (25. 5) 
1 +2acos(w) +a 
(-a)* 

Ry (2) ae R > 0 (226) 
1 -a 


Appendix C shows the results of these equations 
graphically for various values of a. The relationship be- 
tween the power spectra and their corresponding autocorrela- 
tion functions conforms to the expected results from theory 
(l.e., low frequency spectrum with smooth autocorrelation 
function, and high frequency power spectrum with rapidly 
varying autocorrelation function) [Ref. 6:pp. 139-142]. The 
plots in Appendix C also demonstrate that: 1) For poles on 
the positive side of the real axis in the Z plane low fre- 
quencies predominate and for poles on the negative side of 


the real axis high frequencies predominate, and 2) Lower 


a 


magnitudes of a result in a broader power spectrum and a 
wider range of frequencies of significant magnitude. Both 
of these observations agree with the image results. For 
images generated using a more negative pole in a given 
direction, fine, high frequency graininess is observed in 
that direction (though the low contrast or low variation 
about the mean intensity may tend to make this effect less 
noticeable). When a more positive pole is used, lower fre- 
quency variations are more evident in the corresponding 
direction. As lower magnitudes are used for a in a given 
direction, more random variations (indicative of a wider 
range of significant frequency components) are observed in 
that di rection: The form of the autocorrelation funecuwoen 
for these cases approaches the autocorrelation function for 
white noise, 1.@€., an impulse. Negative poles should yield 
high frequencies since the negative side of the real axis in 
the Z plane represents a digital spatial frequency of1, 
while positive poles in the Z plane correspond to a digital 
Spatial frequency of zero. Note that even when the poles 
are placed such that a spatial frequency of zero should pre- 
dominate, there are some low frequency random variations in 
the resulting images. Since the power spectra of the posi- 
tive poles all contain some non-zero frequency components 
(they are not perfect impulses at zero), this characteristic 


1s expected. 


ZG 


C. FOUR POLE, SEPARABLE AUTOREGRESSIVE MODEL 

For the cases considered in this section, H(z) again can 
be factored into expressions in z, and 25. However here 
each factor is a 2nd degree polynomial with two poles in the 
denominator. Figure 2-10 illustrates the filter structure, 


the applicable difference equation, and the corresponding Z 


transform. 
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y(n,m) = - ) ) a, .-y(n-1,m-j) +w(n,m) 
i=0 j=0 ~J 
(1,3)4 (0,0) 
miz.,z2,.) = H(z z.) -W(z a a ee we Ze) 
ee jv? oe l+a 2 sta i l+a gmt+a 272 pre 
HOSE 20st Oia OZ=2 
BE. ee ne ee cy (z Ze) 
~2, za aaa aN -1, ~2, -2 -2 ea 


Ira, 92 #25924 Leg -> 9 ty o2 Pee ee 90-12 


where 


PTO mo a 19 “O22: 712 ~ 702° “10' 722 7 %20°%o2 


Figure 2-10 Autoregressive Filter Form, Difference 
Equation, and Z Transform 


Ca 


Since all of the aj4 coefficients are real, the poles 
must 1) both be on the real axis, or 2) occur in complex 
conjugate pairs in the zZ2, and Zp planes. Again, pole magni- 
tudes must be less than 1 to ensure filter stability. We 
will assume here that the poles in each of the factors have 
equal magnitudes and opposite (or 0 or +7) phase. Letting 

41 = Magnitude of poles in the 2, plane 


8, = pole angle (phase) in the 2, plane 


a9 magnitude of poles in the Z5 plane 
95 = pole angle (phase) an tChew2. senane 
and using Euler's relation, the denominators Of H(Zajem 


H(Z5) can be expressed as follows; 


36 -j6 z a 
-] —2 _ 1 -l i ee pL Zee 
l+a) 121 +a5 921 = (1-0, € Z5 ) (1-a,e Z4 ) =] 2c, cos (84) 24 +1424 
36 mae 2 z 
—l —2 _ 2 -l 2a l 2 
1+a9325 tay 525 = (1-0.,€ 2.5 ) (1-a,e Z5 }=al1 20.,008 (8.4) Z. ta5Z.5 
Hence: 
A190 = 72%1,°cos (9 7) ag, = —242°cos( 8) 
20 alr 02 Z 


This gives a relationship between the pole magnitude and 
angle in the Z domain and the filter coefficients in the 


spatial domain. 
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1. Complex Conjugate Poles 


Figures 2-11 #£=‘through 2 =e illustrate images 
generated using this model for various complex conjugate 
pole combinations in the 2, and 25 directions. In comparing 
each of these image textures in terms of the relative effect 
of pole positioning in each direction, the following obser- 
vations can be made: 


1) For images generated using poles of equal magnitude 
and angle in both directions, graininess with no 
directionality to the pattern resulted. Higher pole 
angles yielded finer (higher frequency) graininess and 
less contrast. Lower magnitude poles yielded a more 
random and less structured graininess pattern at the 
same pole angle. 


2) Using a pole angle of zero (pole on positive real 
axis) in one direction and a pole of some non-zero 
angle in the other direction yielded images with 
highly directional sinusoidal patterns. The direction 
corresponding to the pole on the real axis was not 
totally devoid of variation, but variations were slow, 
1.e., of very low frequency. The spatial frequency of 
the sinusoidal pattern can be increased by increasing 
the pole angle. Large magnitude, high pole angle 
combinations yielded much cleaner and more structured 
textures than low magnitude, low pole angle combina- 
tions. Lower magnitude, high pole angle combinations 
yielded less structured textures where directionality 
was evident but the sinusoidal pattern was obscured. 
Low magnitude, low pole angle combinations yielded 
very random, non-structured textures of relatively 
high contrast. 


3) Using poles in the 2, and Zp planes with the same 
magnitude but different pole angles resulted in some 
directionality if there was a sufficiently large 
magnitude and difference in the pole angles. As 
observed earlier, the direction with the higher pole 
angle had the higher spatial frequency. Large pole 
magnitudes (close to 1) resulted in more structured 
but rather low-contrast images (the low contrast 
seemed to obscure the high frequency nature of the 


pattern somewhat). Pole angles in 2, and Z> that were 
of close value made it difficult to detect the higher 
frequency (higher pole angle) aarece ron. Low 
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Figure 2-12 Images Generated Using a Four Pole 
Autoregressive Model (Poles Listed 
Below Image) 
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magnitudes basically negated the pole angle effects 
and yielded a very random, unstructured, high contrast 
texture. 

4) Where pole magnitudes were close and pole angles were 
different in the 2; and Z> planes, some directionality 
in the texture was observed. Again, the high pole 
angle direction yielded the highest frequency. When 
the pole angles in both Z domains had similar values 
and the magnitudes of the poles differed, graininess 
with little or no discernible directionality or 
structure resulted. 

These observations are consistent with classical 
pole-zero frequency response analysis [Ref. 7:pp. 323-331]. 
There is a direct relationship between pole angle and 
spatial frequency ina given direction, and the magnitude of 
the poles affects the amount of structure and definition of 
the sinusoidal pattern of a given frequency in a given 
direction. Higher magnitude poles result in a narrower 
bandwidth of the filter and yield more structure and 
Sinusoidal pattern definition. Low pole magnitudes give the 
filter wider bandwidth and yield images with less structure 
and definition and more randomness in a given direction. 
While directional dependencies are evident given pole 
magnitude and angle in a given direction, it does not appear 
that a pattern in one direction is totally independent of a 
pattern in the other direction. This would be expected, 
even though the model is separable, due to the cross terms 
in the filter structure. 

Filter power spectrum and autocorrelation analysis 


can be conducted in this case, as in the case of the two 


pole model. The derivations for Sy () and Ry (2) are 


aS 


somewhat more involved, and are given in Appendix D. The 


resulting expressions for Sy(w) and Ry ( %) from Appendix D 


are: 
S(O) : (257 
U) OT ana . 
2 1401842 (02-2 [2c] cos (6) cos (w) +a” (cos (26)+c0s (2) ) ) 
0” _cos(26) cos ( (242) 8)-a“cos (26) (2 > 0) 
RVC) a ee a ca (2.8) 
2sin 6 1-a 1l+a —-2a cos (20) 

The plots of these functions for the various pole 
magnitudes and angles used are given in Appendix E. The 
6 = 0 case is equivalent to having 2 poles on the real axis 
at a given magnitude in the Z plane. As would be expected, 


the power spectrum for each model showed higher magnitudes 
at digital frequencies close to the pole angle. Higher pole 
magnitudes yielded sharper, more well-defined power spectrum 
magnitudes at the given frequency, and lower pole magnitudes 
yielded less well defined more spread-out power spectra. 
Low pole magnitudes almost completely obliterated evidence 
of low frequency power spectrum components, and degraded its 
definition and sharpness at higher frequencies. This 
corresponds to the observed results in the image textures 
generated. The autocorrelation functions also reflect the 
appropriate relationship to the power spectra as outlined in 


the discussion for the two pole case, i.e., greater 
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variation in the autocorrelation function indicates greater 
variation between pixels a given distance apart, which in 
turn implies higher spatial frequencies. 
2. Two Real Poles 

Rather than using complex conjugate pole locations 
to obtain real filter coefficient values, two poles on the 
real axis may also be used for a given direction. They may 
be placed at different locations on the real axis, or they 
may be placed together. The latter situation is equivalent 
to the 6 = O (or 6 = +7 if placed on the negative real axis) 
case, as mentioned above. For the two real pole case, the 


relevant equations are: 


- Ih 1 
eae See =] = ==. —-|. <2 
(1-az ) (1-2 ) 1- (ato, )z FOL OZ 
_ alt 
* -l -2 
lta),z ta y5Z 
where 
Ago. = (aa + Op) 
ELI pee acl 


For these experiments a transfer function of the 
complex conjugate pole form was used for the z, direction, 
With o, = 0.9 and 6, = 0. For the Z5 direction a transfer 


function with two poles on the real axis was used. The 


30 


image textures that result for various values of ® 2 anda 
are given in Figures 2-14 and 2-15. Some observations can 
be made about these results: 


1) With poles placed at the same value on Che Z>.reau 
axis, rather unstructured, low frequency variations 
are observed in the Z5 direction. The more positive 
poles result in very slow variation in the image 
texture, while the lower magnitude positive poles show 
more variation in the Z5 direction. When the poles 
are placed in the same location on the negative side 
of the real axis, a low contrast image with some 
noticeable high frequency variations results. 


2) As the poles are moved farther apart on the )Z55geea 
axis, high frequency variations with increasing struc- 
ture and oscillatory form are evidenced in /Ehemma, 
di Eectronm. 

3) When poles with equal magnitude and opposite sign are 
used, fairly structured high frequency variations are 
evidenced in certain areas of the image, while low 
frequency variations are evident in other areas in the 
25 OREe cron. Higher magnitude poles yield more 
discernible, structured variations, while lower 
magnitude poles of opposite sign yield discernible but 
non-oscillatory high frequency variations in certain 
areas of the image. 

Of particular interest is the fact that two poles 
placed at the same value on the negative real axis in thomae 
plane yielded some high frequency variations. This is in 
keeping with the fact that values on the negative real axis 
correspond to a pole angle (and corresponding digital 
frequency) of @=T7. The presence of poles on the negative 
Side of the real axis of the Z5 plane seems to give rise to 
the high frequency variations with gradually more structure 


and oscillatory appearance as the pole is moved to the left 


(more negative). 
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Figure 2-14 Images Generated Using a Four Pole 
(Two Real Poles) Autoregressive Model 
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Figure 2-15 Images Generated Using a Four Pole 
(Two Real Poles) Autoregressive Model 


38 


The expressions for the power spectrum and autocor- 
relation function for the random process produced by driving 
the filter of Figure 2-10 with white noise are derived in 
Appendix F. It is shown there that the power spectral 


density and the autocorrelation function are given by: 





1-2 (a +o, toro, torch) Cos (W) +2004 Cos (20) +0, 
IL 
+ 201, +04 OL 04, (249) 
2+2 241 Q+1 24+2 
oO O Te 
Ry(%) = 1 na _ fa “bvalb  ,b ] (21,0) 


2 Z 1-a Z 
(A,-O) 1-a, ab 1-0, 


Plots of these functions for the various values used 
in this section are given in Appendix G. The power spectrum 
results are consistent with the observed image spatial fre- 
quency characteristics. Both low and high frequency 
components were contained in some of the power spectra, and 
were manifested in the corresponding images as both low and 
high frequency variations in the Z5 direction. The nature 
of the autocorrelation functions related to the power 
spectra that contained low and high frequency components was 
interesting. Autocorrelation functions with much variation 


but all positive values, rather than the equal magnitude 
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positive and negative values evidenced in earlier results, 
seems to reflect the higher level of correlation related to 
the low frequency (smoother variations) aspect of the image 


texture variations. 


DD: IMAGE TEXTURE ROTATION TRANSFORMATION 

If an image signal x(nj,N9) consists of a _ rotated 
version of another image w(m,,M 7) such that n, = Im) + me 
and n> = Km, + Lm», where I, J, K, and L are integers and 
IL-KJ # 0, then the 2 transform X(2Z1,25) is) Sai Venue 
w(zt,z5, 24,25) (Ref. 3:p. 182]. A 45° rotation corresponds 
to Il = 1, ea eee 1, L= =i. If we use the four pole 
separable result for H(21,Z 2), as shown in Figure 2-10, and 
apply the above transformation (2, > 2, Zoe = 2125"), we 


Lind atter Simpl rE1cacion: 


1 
Hp (2,23) = Sp. SCO eee 


=a 7) i 7) 
144,12) 25 +8992] 29481921 22 t210%0171 7210%0271 72 


i" 2.992] 29 *,920121 22 *20%021_ oe 

Notice that this transfer function is not separable but 
consists of a rotated version of a separable filter. Figure 
2-16 illustrates the support of the denominator polynomial 
for this filter. It has the form of a non-symmetric half- 
plane infinite impulse response (IIR) filter, so it is 


recursively computable. 
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Figure 2-16 Rotation Transformation Filter Form 


The application of this filter, using filter coefficients of 
oe 

the four pole separable filter with poles at Zz, > 0.9e 2. 

25 + 0.9e7J9 in the Original separable filter yielded the 


result shown in Figure 2-17. 





Figure 2-17 Result of Rotation Transformation 


E. SUMMARY 
Autoregressive models can produce a variety of image 
textures. For general two-dimensional models, the system 


functions are generally not factorable and singularities 


4l 


occur on surfaces, not at isolated points. For these 
reasons it is difficult to design two-dimensional filters 
for images and predict the resulting character of the 
images. Indeed, even to ensure stability of the filter is 
not trivial. As a result we concentrated here on separable 
forms, which by their nature are much easier to analyze. 
Certain types of texture patterns uSing various separable 
autoregressive models can be predicted based on filter pole 
placement in the Z, and 25 planes. Arbitrary or random 
selection of filter coefficients can yield interesting but 
generally unpredictable results. Obviously, an infinite 
number of variations on the models above could be attempted. 
Ultimately, the anticipated utility of the textures 
generated will guide the process of model and parameter 


selection. 
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ale. IMPLEMENTATION OF AN FIR SUMMATION FILTER 
IN TWO DIMENSIONS 
To implement the 2-D ARIMA model, the inverse of a 
filter representing a suitable difference operator is 
needed. One possible 2-D difference operator is the 
Laplacian, which has the impulse response shown in Figure 3- 


MeiRet. 8:pp. 212-213]. 


m ney, mM) 





Figure 3-1 Laplacian Impulse Response 


Its implementation involves convolving it with an image 
and is represented by the following difference equation: 
1 dl 


y(n,m) = ) ee ae eal) (38) 
i=-1 j=-1 


where 
x(n,m) is the image input signal, and 


fowler am Gocrh Clent matrix 


b:- 
“7 (ba5 = (G3) 
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In the Z-transform domain this can be written as: 


¥(24,29) = BZ] 7 20 eee 


(4 - z71 - Zo a Z9) °X(Zy;7 Zon 


In areas of an image where adjacent pixels have similar gray 
levels (low frequency, homogeneous areas), the result of 
this operator will be approximately zero. Where significant 
or sharp differences in gray levels between adjacent pixels 
exist, the result of this operation will be farther from 
zero. Thus the Laplacian difference operator is sometimes 
used as an "edge detector." 

The problem addressed in this chapter is constructing 
the inverse of the operator. In the Z domain, the 
expression for the inverse would be [Ref. 4:p. 36]: 

es i i 
which has an expansion as an infinite series of positive and 
negative powers of Z, and Z>. That is, considering Jem 
expression as a problem in long division, the result of such 


division would be an expression of the form: 
H-1(z4,29) = A ) a Sel es 


where 
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= coefficient values of ztz4 resulting from the 


4 ee 
long division 


Note “ehateirf 





1 * co CO epee 
=I -1_5. =a ib aa 
4 - Z, = Z4 Zo 5 J 
then cross multiplication would yield: 
eb z—1 2) ae z~2 = Ze , ) ee ecae = 1 (37.3) 


The double summation expression in 2, and 2Z 9 will be 
truncated and considered to be an FIR filter with finite 
support and coefficients ais: This approximates the desired 
inverse filter. In particular, we will use the following 
constraints: 

1) Choose the limits of summation to be equal in both 


directions, il.e., 


th =o “i 
te ss rile 


This results in a "square" region of support for the 
filter (all values outside assumed zero). 

2) Force the values for the filter coefficients to be 
Symmetric, 1.e., ajy = a-jj = Aaj-j = Aa-j-j = A4i 
mde > Aq=q = aay -4- 

Using these constraints and implementing the cross 
multiplication equation (3.3) will result in an expression 


in Z, and z5, with each combination of the zizJ, terms having 


a coefficient whose form is a summation of terms in aj; 
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where the coefficients of aj4 are either 4 or -1. The 
coefficient of the 2920 term must equal 1 and the coeffi- 
cient of any other Ze term must equal zero to satisfy 
equation: (3.30). 


As an example, let L= 1. Equation (3.3) can theneee 


expressed as: 


il 1 
(4 - 2-1 - 2, - 271 - 25): a wat 21571 om 


2 Jlheeak 
= jt 0°24 + 0°25 + 0325/25 eee 


Performing the double summation yields 


i il ie pee 
»-sg7loam s lol 1,0 1o-1 
| a )  a452] ae = €-7-127 23 + 48-1027 227 + a-112725° + 
i==—1 j—— : 
0,0 O71 -l,l =i 
+ Ag92°2~ + ag- 12-2, + aj—y42g—2, tee 
ors 0-12, 2, 1-12) ~25 1021 


+ age 7 442) 7257 
Performing the cross multiplication would yield 45 different 
terms in various combinations of Zee Combining these 
terms to find the coefficient expression for each Za aoe 
term soon becomes rather tedious and impractical for even 
moderate values of L. An alternative way to proceed is to 
choose a A term on the right hand side of the equation, 
and for each term in the expression sa25b=24 cee eee deter- 
mine what values of i and j are required so that when 
each term is multiplied by ayyz5i2zy. it will result in an 


Z 
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expression in the chosen 2228 term on the right hand side of 


the equation. Choosing 202. (= 1) on the right hand side of 


the equation, we have: 


Bn sez 5S = ez when i = 0 and j = 0; so c = 4aQ0 
-27t-ajyz7425) = C2020 when i = -1 and j = 0; So c = -Aa.10 
B24 A442, 25° = C2025 when i = 1 and j = 0; so c = -ajo 
eee aq52, 725) = cz 0 20 when 1 = 0 and j =-1; so c = -aoq 
“29° a442, 425) = cz0 20 when i = 0 and j = 1; so c = -ao) 


So the coefficient for 2 1s simply a summation of the c 


terms obtained above, i.e., 


= = = 050 
(4agng9~ 410 ~ 810 ~ ao-1 ~ 801) 2525 
This entire expression must equal 1 to satisfy (3.3), and 
since Z2 20 Tee fea) = ayy Ag=] — ag, — Ll also. 
Using the same method for the 272. term on the right 


hand side yields: 


4°ai52, z~J = CZ 72% when 1 =-1 and j =-1; so c = 4a_}-}] 
ae meciyes 2) = C2521 when mH Eaceande je —l;, SO C — —al>_7 
21 aunZe J= cz>27, when 1 = 0 and j =-1; so Cc = -ap-} 
Z Yeaqj27°27 = cz721, when aioe —2; co C — -a_7 > 
29 aaa = cz727, when Peo leond) 3 — "0; SO C = —al7, 
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The resulting term in Zea is: 


(4a-y-1 ~ 4-2-1 ~ 89-1 ~ 4-1-2 7 a-19) 2425 


This expression must equal zero, since there is no z 2 term 


af 
i 
on the right side of (3.3), so 4€a_,-4 = €@=52]q) = aq =n 
== Qn ke 

This procedure can be extended to any number of 25 ee 
terms. When this is done, the resulting expressions for the 
coefficients of 274273 can be formed into a set of simultan- 
eous equations in order to solve for the aj4j coefficient 
values. However, due to the symmetry condition imposed 
above, some of the equations for the coefficients of the 
Zee terms are linearly dependent. For values of i and 3 
that yield unique or distinct values for aj4, the resulting 
Zee coefficient expressions are linearly independent. 


1 
For example, the coefficient expression for the z—1z92 term 


12 
is linearly independent of the coefficient expression for 
the cae term, since a10 * 411- But the coefficient ex- 
pression for the 252 term is linearly dependent on the 
coefficient expression for the Ze term, since aj19 = 491- 
Using only the linearly independent equations for a given 
filter size yields a set of p equations in p unknowns, where 
p is the number of unique and distinct ajj values in the in- 


verse filter. The value of p is related to the size of the 


desired inverse filter. If the size of the filter is Nias 
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the number of unique aj; values using the symmetry of 
constraint above is: 
aNd 


A +1) + = + SS -1) eal (3.4a) 


se. 
l 


= ——_— (3.4b) 


For example, with N= 7, the unique aj4j values in a 7x7 
iaverse filter can be represented by 433, 839, 431, 430, 
elo, o>o7 4117 4107 200. “hough there are 49 elements 
in a7x7 filter, all of them are equal to one of these 
values listed, due to symmetry. Obviously, N is constrained 
to be odd for a square filter with a unique element ago in 
the center. 

The solution of the resulting p equations yields the 
values for the p filter elements or coefficients. This 
defines the FIR approximation to the inverse filter. It is 
only an approximation due to the finite size constraint 
imposed, and it might be expected that the larger the filter 
Size, the better the approximation. 

An algorithmic procedure for obtaining the ai5 
coefficients is outlined below. An example follows. 

1) Determine the desired size of the inverse filter. 
2) For each combination of (positive) i,j values corres- 
ponding to a unique ai; filter coefficient, identify 


the five term summation equation associated with each 
ZqjtzzJ term. 
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3) Combine equal aj;+ values and develop a matrix of 
coefficients for the ajj values. Let this matrix be 
A. 


4) Denoting the vector of unique a;~ values as a, the set 
of simultaneous equations in matrix form is: 


| 
So 


Aart = (3.5) 


LJ 


where 


a = [890 410 411 4820 421 4422 -:: 2y-1 wed 


The top row of A corresponds to the summation of ais 
terms that represents the coefficient of the 2929 
IL 
term. 
5) Solve (QS) tfor ae 
An example of this procedure iS appropriate at this 
point. For an inverse filter of size NxN: 
Step 1 
Let N = 5 (therefore L = 2) 
step 2 
The coefficients corresponding to each unique z-iz7) 


Term awex 

5 «6 «= 4800 - 810 - a-10 - 401 - a0-2 
25 * 4810 — 820 ~ 400 0= (2 ihieieyacnees 
1 


17% > 4411 7 421 > 4OnMie 4r2ieeno 
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aa OTe Ou nOme cD] ~ A22—7 


Ze eo tees ewe — 229 — ap—-2 
i i ee OG ag) Dea 3 AD 
Step 3 


Combining equal terms in Step 2 and expressing the 


coefficients of aj4 in matrix form yields an A matrix of: 


4 —-4 O 0 0 0 

-1 4 -2 -1 0 0 

—— O -2 4 0 -2 O 
Oo -1 0 4 -2 0 

0 O -1l -1 4 -1 

0 0 O 0 -2 4 


Step 4 


With A given in Step 3, the a vector for (3.5) is 


€a = [80 410 411 4820 421 422] 


Step 5 
Solving (3.5) for a involves inverting A and multiplying 


it by [1 0 0 0 0 ojf. Thus: 


Dl 


Appendix H illustrates the forms of the resulting 
inverse filters of various sizes, as well as the normalized 
and unnormalized filter cross sections. 

One way to validate the resulting inverse filter is to 
convolve it with the original Laplacian difference operator. 
The result should approximate an impulse at the origin. 
Appendix I shows the results of this convolution using 3 x3, 
5x5, 7x7, 9X9, 15%* 15, and 21x 21 size inverse filters. 
It is seen there that as the size of the filter gets larger, 
1t becomes a better approximation to the true inverse and 
the convolution looks more like an impulse. 


To test the application of this filter on an actual 


image, a test image was filtered using the Laplacian 
difference operator. Then the resulting image data were 
filtered again using various size inverse filters. The 


results are shown in Figures 3-2 and 3-3. Note that the 
image resulting from Laplacian FIR filtering seems more 
stationary than the test image, which was one of the desired 
results. Inverse~.filtering of that result yields images 
that are progressively more similar to the original test 
image as the size of the inverse filter increases. However, 
a rather large inverse filter is needed to accurately 
reproduce the image. The result of the 21x21 inverse 
filter is quite similar to the test image, with some loss of 
contrast or darkness in certain areas, but with essentially 


the same pattern. The effect of the size limitation of the 


DZ 


. CRS ROE LSA SIS Re 
: SY RO NPS OS OAS 





Test image Laplacian filtered test 
image 





3x 3 inverse filter 5x5 inverse filter 


Figure 3-2 Results of Filtering Test Image with 
Laplacian FIR Filter and its Inverse 
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Figure 3-3 Results of Filtering Test Image with 
Laplacian FIR Filter and its Inverse 
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inverse filter, as manifested in the convolution of the 
Laplacian and its inverse above, would seem to explain the 
lack of perfect test image reproduction. Larger inverse 
filter sizes could be tried, but large inverse filter sizes 
relative to image size would result ina significant portion 
around the edge of the image having only a part of the 
filter applied to it. This would adversely affect the 


overall quality of image reproduction. 
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IV. APPLICATION OF THE ARIMA MODEL TO IMAGE TEXTURES 


As outlined in Chapter I, the utility of the ARIMA model 
centers around the fact that a difference operator applied 
to an image texture may improve the stationarity of the 
image statistical characteristics. A stationary image 
texture is required for accurate modeling by autoregressive 
techniques, and it was hoped that application of the 
autoregressively generated texture to an approximate inverse 
of the difference operator may yield a more accurate or 
recognizable representation of the original nonstationary 
image, aS compared to a purely autoregressively generated 
version. 

A. APPLICATION OF LAPLACIAN INVERSE FILTER TO AUTO- 

REGRESSIVELY GENERATED IMAGES 

As an initial examination of the effects of the inverse 
filter developed in Chapter III on image textures, selected 
images generated in Chapter II were input to the 21% 21 
version of that filter. Figures 4-1 through 4-3 illustrate 
the results. All attempts resulted in a blurred or smoothed 
version of the original image. Since the inverse of a 
difference operation is a summation or "integration" 
operation, and since integration operations can be expected 
to smooth or blur (low pass filter) signals [Ref. 8:pp. 136- 


LS Aa the results are Met surprising. However, 
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Figure 4-1 Images from Figure 2-7 (top 4) Applied to 
Summation Filter 


a7 








Figure 4-2 Images from Figure 2-11 (top 4) Applied to 
Summation Fiver 
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Figure 4-3 Images from Figure 2-14 (top 4) Applied to 
Summation Filter 
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except to the extent that blurring is useful or desirable, 
applying the summation filter to image signals that are not 
based on the application of the corresponding difference 
filter to that signal seems to be of little utility. 

In the remainder of this chapter we consider applica- 
tion of the summation filter to regenerate actual image 


textures. 


B. AUTOREGRESSIVE FILTER PARAMETER ESTIMATION PROCEDURES 
The first step in testing the ARIMA model is to estimate 
the autoregressive, quarter plane filter parameters required 
to model the real image textures and the signal resulting 
from application of the Laplacian operator to those images. 
For a zero-mean signal, these model parameters are found by 
solving a set of Normal equations. In these equations the 


white noise covariance is referred to as the prediction 


error covariance. The Normal equations can be expressed as 
a = 
Me ay 0 

R ® = ® { 4: at) 
ao= 1 0 


where the R matrix is the correlation matrix for the signal 
(block Toeplitz with Toeplitz blocks), the a vector consists 
of appropriately ordered filter coefficient vectors, and s 


is a vector containing the prediction error covariance as 
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the first and only nonzero element. Here aj = hang Aq] 
Teco). and S$ = [6° 00. . #O]t 

Calculating the correlation matrix and prediction error 
covariance from the image signals, and solving (4.1) for the 
a vectors, provides all the parameters needed for the 2-D AR 
model. The multichannel form of the Levinson recursion can 
be used to solve these equations more efficiently [Ref. 2:p. 


454). 


C. APPLICATION TO REAL IMAGE TEXTURES 

Actual image textures used here are from the image data 
base at the University of Southern California's Signal and 
Image Processing Institute [Ref. 9:pp. 13-14]. The images 
selected from this data base are contained in a book by 
Brodatz [Ref. 10]. Portions of the images of size 128 x 128 
pixels were obtained and used as a basis for processing. 
Filter coefficients were calculated for the real image 
textures shown in Figures 4-4 and 4-5. Various filter sizes 
were tried to determine which yielded the best results in 
generating a particular image, and a quarter-plane filter 
size of 4x4 was’ selected. Results of autoregressive 
filtering of white noise using the appropriate calculated 
coefficients to model each texture are given in Figures 4-6 
and 4-7. Images generated by applying the Laplacian differ- 
ence operator to the real images are shown in Figures 4-8 


and 4-9. Autoregressive generation of these images using 
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Figure 4-4 
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Figure 4-7 Image Textures Generated Using an AR Model 
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Figure 4-8 Actual Images After Laplacian Filtering 
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4~9 


Figure 


filters with the corresponding calculated coefficients are 
given in Figures 4-10 and 4-11. Finally, the application of 
the signal represented by the images in Figures 4-10 and 4- 
11 (without the 0-255 scaling reflected in these figures) to 
the 21x 21 inverse filter described in Chapter III yields 
the images shown in Figures 4-12 and 4-13. Comparison of 
all of the above results yields the following observations: 

1) Autoregressive modeling of the water, grass and sand 
textures yielded good results. Some of the other tex- 
tures with more structure and sharp local variations 
were not reproduced well. 

2) Autoregressive reproduction of images created after 
application of the difference operator, with the 
exception of the water image, yielded generally poor 
results. As observed in Chapter III, the application 
of the difference operator produces a seemingly more 
stationary result, but the edge structure that remain- 
ed in most of the images after application of the dif- 
ference operator was in general not reproducible using 
a purely AR model. 

3) Application of the inverse filter to the image signal 
generated by AR model reproduction of the difference 
operator results yielded smoothed versions of those 
results. This is similar to what was observed in Sec- 
tion A of this chapter when images were applied to the 
inverse filter that were not based on the specific 
data generated by the difference operator. 

As a final test of the ARIMA model, a 64 x64 contrast 
enhanced aerial photograph oof trees, with smoother 
variations and in general less edge structure than the other 
images tested, was tried. The results are shown in Figure 


4-14. Though this image seemed somewhat better adapted to 


the model, overall the same observations outlined above 


apply. 
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Generated Using an AR Model 
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Figure 4-11 Laplacian Filtered Image Textures 
Generated Using an AR Model 
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Figure 4-14 Final Test of ARIMA Model on Contrast 
Enhanced Trees (Magnification X2) 
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De SUMMARY 

The effectiveness of AR reproduction of image data using 
white noise input and filter coefficients calculated based 
on the statistics of the image signal is highly dependent on 
the nature of the image data. The water image, for example, 
with its smoothly varying and rather homogeneous nature, was 
quite well adapted to AR reproduction. Other images with 
more structure, abrupt variations, and more non-homogeneous 
characteristics, were not autoregressively reproducible to 
any great extent. 

Using the ARIMA model, it seems that the operation of 
the inverse filter is very sensitive to the nature of the 
input data. Input data that are strictly based on the 
difference operator output can reproduce the original image, 
as was found in Chapter III. However, the AR model used to 
generate the inverse filter input (based on the statistics 
of image signal produced using the difference operator) does 
not generate image data accurately enough to reproduce 


images that resemble the real images tested. 
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V. CONCLUSIONS 


This thesis sought to explore experimentally and to 
understand how linear filtering models could be used to 
generate texture in images. Of particular interest was the 
investigation of 2-D ARIMA models to see if they might be of 
any utility in this effort. Some time was spent exploring 
separable 2-D models to understand how transfer function 
pole placement affected image texture characteristics. 
Image textures generated using these models and applied to 
the summation filter yielded blurred or smooth textures with 
seemingly little variety or utility. The ultimate test of 
the model was its ability to reproduce actual image 
textures. The purely AR portion of the model reproduced a 
few types of actual textures well. However, the full ARIMA 
model failed to generate image textures that resembled the 
source images used. Many of the textures had strong edge 
differences that were not accurately reproducible by the AR 
model. Also, the summation filter developed seemed very 
sensitive to deviations in image signal data from that 
generated by application of the difference operator; that 
is, the procedure seemed not to be "robust." 

Since many of the images tested here have definite edge 
structure, the difference image had lines which were not 


reproduced well by the AR model. Correspondingly the 


uo 


integrated AR model did not reproduce the original image. 
For this type of image, a combination of a line point 
process model [{Ref. 11] with the integrator, would possibly 
have been more suitable. The image of trees had not such 
edge structure and produced somewhat better results. 
Further experimentation with images of this type and the 


ARIMA model would perhaps be appropriate. 
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APPENDIX A 


} 
COMPUTER PROGRAMS, SUBROUTINES, AND FUNCTIONS 


Listed below are the names, associated computer systems, 


and functions of the various computer algorithms used to 


accomplish this thesis research. All programs were written 


by the thesis author unless otherwise noted. Program source 


codes are given at the end of this appendix (except for the 


MAKFIL* series). 


A. 


1) 


2) 


3) 


4) 


>) 


PROGRAMS 


AUTOREG (VAX/VMS FORTRAN) 


The program did the following: 


Generated a 128 x128 zero mean white noise matrix 
using subroutine PGAUSS. 


Multiplied the white noise by the appropriate image 
data standard deviation when necessary. 


Converted that matrix into a displayable image file 
using subroutines SCALE and INTBYTE, when necessary. 


Read filter parameters from an input file into an 
array. 


Implemented the equation: 


Ose 
y(n,m) = - Y “Y ayyy(n-i,m-j) + w(n,m) (A-1) 
1 


using the white noise array and the filter coefficient 
array as inputs. 


ey 


6) Converted the array result from 4) into a displayable 
image file using subroutines SCALE and INTBYTE. 


7) Used subroutine SUBINTFILE to create image data files 
from filter results for further processing. 


8) Used subroutine NONC to apply the summation filter to 
image data when necessary. 


2. NONCAUSAL (VAX/VMS FORTRAN) 
This program did the following: 
1) Read filter coefficient values into an array. 
2) Read image data from an input image file, converted it 
to integer values using subroutine BYTEINT, calculated 
the mean from the data, and placed the data into a 


real array. 
3) Implemented the equation: 
L L 


y(n,m) = y y aj4x(n-1i,m-j) (AS 28 
=o 


using the image data array and the filter coefficient 
array. 


4) Called the subroutine NONC to implement the equation 
in Step 3) a second time, when necessary. 


5) Converted the result of Step 3) to a displayable image 
file using subroutine SCALE and INTBYTE. 


Se CONV _(VAX/VMX_FORTRAN ) 
This program performs the same basic functions as 
NONCAUSAL, without having the capability of calling 
subroutine NONC. It was used for convenience in convolving 


certain filter structures with certain test images directly. 
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4. MAKFIL* (VAX/VMS_ FORTRAN) 


This family of programs was used to create various 
autoregressive and FIR filter coefficient files, using 
source data manually entered into the program. 

Dee oe ECORZ (TBM SYSTEM/3m0n 3033. VS FORTRAN 1.4.1 

This program implemented the equations derived in 
Appendix B and created data files used in developing the 
corresponding graphs. 

Gn OPE CORs eI BMeSYSTEM7 37023033. VS FORTRAN 1.4.1) 

This program implemented the equations derived in 
Appendix D and created data files used in developing the 
corresponding graphs. 

Rice SEepCOkRSA( LUBMvsYoTEM/370 3033 VS FORTRAN 1.4.1 

This program implemented the equations derived in 
Appendix F and created data files used in developing 
corresponding graphs. 

8. VARIMGS (VAX/VMS FORTRAN) 

This program was used to display image data files on 

the COMTAL (not written by author). 
9. PIECE (VAX/VMS FORTRAN) 

This program was used to make 128x 128 image data 

files from larger image data files. 
10. INTFILE (VAX/VMS FORTRAN) 
This program created appropriately formatted integer 


files from input image data for further processing. 


a 


11. TRANS (VAX/VMS FORTRAN) 


This program changed the format of filter coeffi- 
cient data files into a form readable by the image 
processing programs. 

12. NSHP (VAX/VMS FORTRAN) 

This program was used to convert quarter-plane auto- 
regressive filter coefficient data to non-symmetric half- 
plane autoregressive filter coefficient data based on the 


transformation outlined in Chapter II, Section-D. 


Be SUBROUTINES 
i. PGAUSS (VAX/VMS FORTRAN) 

This subroutine, written by C.W. Therrien, was used 
to generate zero mean, unit variance white noise using RAN 
(a random number generator function) SQRT, COS, and SIN 
FORTRAN functions. 

Pee SCALE (VAX/VMS FORTRAN) 

This subroutine takes an image data array and 
converts it to an integer array with values between an input 
Maximum (MAX) and minimum (MIN) using the following scaling 


formula: 


(A(i,j)-LOW) x (MAX-MIN) 


A.3 
HIGH-LOW — (Aaa 


Gli a) 


A(i,j) is the input image data array, I(i,}) is the output 
integer array, and HIGH and LOW are the high and low values 


of A(i,j), respectively (calculated in this subroutine). 
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This is done to provide appropriate values for image files 
that will be displayed on the COMTAL Vision One/20, since 
the gray scale intensity level of each pixel is represented 
by an 8-bit word. So values possible (in base 10) range 
from 0 (darkest), to 255 (brightest). 

a, INTBYTE AND BYTEINT (VAX/VMS FORTRAN) 

These subroutines are necessary since data in an 
image file are stored in two's complement form. The related 
variable type in FORTRAN for these values is BYTE. ro 
process image data using FORTRAN implementation of the 
appropriate formulas, these byte values must be converted to 
integer (and eventually real using the FLOAT function) form. 
Results of image processing formulas in real form must be 
converted to integer (using the INT function) and then byte 
form to be placed in image data files. INTBYTE converts 
integer type variables to byte type variables using the 
following criterion (I is an integer and B is a byte): 

If I < 127 and I > 0 then B= I 

Reel > 127 and £ < 255 then B = 1-256 
BYTEINT converts byte type variables to integer type 
variables using the following criterion: 

If B > -128 and B < O then I = Bt256 


If B O and B < 127 then I = B 


4. SUBINTFILE (VAX/VMS FORTRAN) 


This subroutine performed the same function as 


IV 


INTFILE, but could be called by a program to operate on 
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processed image data arrays, rather than just image file 
data inputs. 
5. NONC (VAX/VMX FORTRAN) 
This subroutine performs essentially the same 
functions as NONCAUSAL, except that it can be called by a 


program to operate on an image data array. 


C. APL FUNCTIONS 
The APL systems on the IBM System/370 3033 and VAX/UNIX 
were used for matrix manipulations and operations, for 
graphing filter structures and convolution results, and for 
calculating autoregressive filter coefficients from image 
data. All APL functions except MAKMAT were written by C.W. 
Therrien. 
1. MAKMAT (IBM 
This function was used to create the large coeffi- 
cient matrices (A) used in calculating the FIR fue 
coefficients as outlined in Chapter III. 
2a CC22 1 BM 
This function was used to circularly convolve the 
Laplacian FIR filter and its various inverses. Appropriate 
zero-padding of these filters makes the resulting circular 
convolution equivalent to linear convolution (Ref. 2:pp. 70- 
72], which was the desired operation. 
3. GETDATA (VAX/UNIX 
This function is used to transfer image data files 


from a UNIX subdirectory to an APL workspace. 


82 


4. PUTDATA (VAX/UNIX) 


This function is used to transfer filter coefficient 

data files from an APL workspace to a UNIX subdirectory. 
5. MEAN (VAX/UNIX) 

This function is used to calculate the mean of an 

image data file for use in the APL function COVF. 
6. COVF (VAX/UNIX) 

This function is used to calculate terms in the 2-D 

covariance function for use in the APL function CORR. 
7. CORR (VAX/UNIX) 

This function is used to estimate the 2-D covariance 

function of the image data. 
8. MVLEV (VAX/UNIX) 

This function is used in APL function FF2DLEV to 
calculate necessary parameters for the 2-D Levinson 
recursion from the covariance function of the image data. 

9. FF2DLEV (VAX/UNIX) 
This function performs the 2-D Levinson recursion to 


solve for the filter coefficient vector. 
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annannianna 


NOnnnAnenankanaanna 


AUTOREG 


THIS PIDSRAM GENERATES AY TMOGE TET JRE YSING s+HIFTe “O7S=E aS av INPUT 
TO AN 4UTOREGRESSIVE FILTE? sHNSE S49A8OME TERS ARE OBTAINED FROY THE 
FILE FILTOEF. SYARJITINVE PSAUSS 1S YSED TO GENERATE RTHE [INPUT Sara 
VOISE 4N) SURRAUTINES SC4LE ANN INTAYTE AR= YSED 19 PREPARE [MAGE 
JQ4TA AIRAVYS FOR DISPLAY. SIJBSOUTINE SUSINIS ILE IS YSEDO IF AN 
IVIEGE? FILE RESULT TS DESIRED, 4ND SUSRQNUTINE NONC IS USED IF A 
LAPLACTIAW INVERSE FILTERING STEP I5 NEEDED. 


DEFINE VARTABLES 
Dyte af 92t27),bim(02127,0:127) 
Integer Np SeeteFSiZercsiZzecOMlsamisisjriOlsrowscoleinregl(02:127.,021 
'27) 
realt3 vall,valeran(O2127,0: 127 )otsun,inl(OsleT»,Qsi27),coet(039,039 
Sy evaremax, HIN, 1 ROUt(O03127,92127) 
OPEN FILES 
open(unit=l,name=' (rarnnaran.datal an.dat'styve='new' ,access='direct 
‘* recor gsi 72e532,maxre-=129) 
onen(unitse,name=*rarnnaanvinl final%a.tat ',tvoe='new's,access='ai 
‘rect '-recorisi1ze=se,maxrrerzsi1ce%)} 
osen(3,file='rathnann.ratal fcdous.sar‘’,status=‘old') 
DEFINE P4RAaMETERS 
seeysl123%4567 
rsize=127 
Csizesi27 
onl=3 
a3n1=5 
max=229,9350 
nminse7.jan 
CREATE aHtITE NOISE ARRAY 
65 #19 itd,rsize 
do 20 jy=0,esize-! 
polsj;el 
call sqauss(seedevali,vale) 
wrlisjdsvall 
wrlisjol)devale 
20 contisaue 
10 €5ntinte 
SC4&Le ARPAy aN ANDO CANVERT TQ BYTE FIR 
call srate(wrr,inteas, rane nin) 
call tatovre(integ,d1m) 
ARTTE THE wHITE NOISE [MAGE ARRAY TO A FILE 
Jo 30 j20,es12ze 
39 42 jy F0-2512e 
aC jpdsorn(t,§) 
an eortwnue 
arerytel(iiet) Calan) ,n20,ersize) 
39 continge 
SEQC FILLIES PASE*ETESS [ld Sw Senay 
J5 30 F20.90%1 
ds 50 j2=9,.2-1 
peoat(3,55) soef(1,3) 


52 fornar (420.12) 
ov continue 
5° ¢antirue 
Ves TPL Y AAITE NDISE By RE&L PURGE SIA OaR DOE eG. 

varzsarr(caef(N,9)) 

35 70 y2U,ersize 

G5 30 j=Q,csize 
wrlisvjdtarnlr, i devar 

89 CoONtirue 
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anNanNnaAaaanNn 


AUTOREG (CONT. ) 


70 contiraue 
aPOLy wHITE NOTS=E TO THE AJTOREGRESSIVE FILTER 
O5 110 r=%,-rsize 
a5 129 m=0-28172e 
osur=0.000 
92 1350 VeH,om) 
On 14% FENn,av] 
it Gty eas0) sand. j.e0.0)) a5 to 140 
FPAe n=} 
calione-j 
it Crowle .-O).ore(enl.le.9)) ad to 140 
Osum=a(coeflics) sinlrowecol))¢asum 


140 continue 
130 continue 

in€(asnde(-1 Od0tasun) ewn (nee) 
120 continue 


110 esnrinue 
FILTER THE IVAGE DATA SRRAY JSING THE LAOLACTAN INVERSE FILTER 


call roaclinegrernut ) 

SCALE THE VZESULTING IMAGE ARRAY AND CONVERT TO BYTE FORM 
call scalel(ynsuterrterasnaxre min) ; 
call inerovtel(inrearnin) 

ARITE THE SENER ATED IMAGE INTO A FILE 
O05 15) it0,rsize 

O5 169 j=O-zsiZe 
al(jreorm(irj) 
160 contiave 
werte(c'ie*l) Calnden=Oersize) 
1S0 conrirue 

Edo= = TLES 
Close(usit=1 ) 
closeluait=2) 

Close(3) 
stos 
end 


NONCAUSAL 


THIJS PVWDIRA4 GENERATES AN FMAGE TEXT YRE SINS & VOUCAUSAL FIR FILTER 
whISE PATALETERS AQE V4lALyEY FROM A DATA FILE, THE FILTER IS 
SOOL ITE Tu At TMAGE. SIJAQDITINES SCALE, INYIVTE, aid SYTEINT ARE 
JSED TS SQFPav=E PMAGE DATA ARRAYS FOF NISPLAar. SJAROUTINE NONC 
15 USED IF Av INTEIPMENLATE FILTERING STE® 1S DESIRED, 


NESLYE VaRTARLES 
Dyt® se():ie7),hbim(O0s1e7,02127) ; 
Vrateger AerstrercCsize@riesjeoPowrcCol,rrates(Osle7,VzleTjetsizercntfsize, 
e1ndexleintgexscc 
real*3 asus,in(Os127 02127 ),coWetfl-10210,-10:10) ,nsumenean, dunrnl 03] 
27 DSICT)AHigns tower tevlaels128, 13122), 2an, nin 
oP2N FILES 
osen(unitel,names' (rathnaran.inivarx«.dat'styoe=‘'old',access=‘daier 
fect 'srecsorasizezser,naxrrectlea) 
s0eNn(UuraAItsr,nare='[ratnmann.inlarindO%. gat',tvoe='new'saccess=' dire 
tot 'ererordsizer32-nanrec=128) 
osen(3,file='fratnnansn.taralacfilcota.sat sestatus='oly') 
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Cc 


c 
C 


NONCAUSAL (CONT. ) 


DEFINE PARAMETERS 
indexl=~-] 
Vydex2=128% 
rsize=1]>7 
csizeslo7 
fsize=] 
mf#Size=-lefsize 
max=235.,030 
njinz0.040 
READ FILTE2 PARAVETERPS INT) AN ARRAY 
O65 10 isfsizer,fsize 
05 20 jrontsize-fsize 
eeai(3,25) co4“#flisj) 
25 forrnat(320.12) 
20 contisue 
10 continue 
READ [MASE TJ) BE FILTESEDS INTO AN ARRAY AND 
O05 30 120,rsize 
eead(i'it*l) €aln),n=0,127) 
d> $0 jr0,csize 
binl(i,j)aalj) 
ud cCantisaue 
30 contirus 
Call svreintloinreintes) 
CONVERT THE INTESE2 ARQAY INTO A REAL ARRAY 
nsurn=).0d00 
hiqgn=-19000,050 
bow2=10099.0790 
95 30 i2f0,rsize 
d>. 5) p=9,¢s12¢ 
iu(ysj DSP lost Crats3( 1,99) 
~sunzin(i,s)¢rSur 
Vo (hin pain ie Se 6 4) NOS alee aged 
Vt Cp wigs 6 tole) Pees ia 
50 CONriyue 
SO contiryue 
mearansym/(12%.9*#12%.,9) 
arice(s,55) nean,hiranelo« 


CONVERT JOS1 NE Gee 


AND COMPUTE THE MEAN 


SS forrar(’ griaingl irnage Jdatarz nearne',3912.5-' higns ,die.s,° 


aol255) 
SURFIATT THE VEAY FRIA THE OTS 
35.70 Veuy l 27 
V0 aml 
me(iseplsrnliey leomear 
80 conrisaue 
70 ¢conti4Us 
FILTER FHE [4AGE array 
Jn 90 nzindenlsindex2 
Q> 104 xvntindencl,indtex? 
3Ssun=J.990 
3> 110 tiznfsize-fsize 
45 129 jenfsize,fsize 
r{yasna-i 
colar-j 


i¢€Cese.lt.9).0r. (col clit.) aaet ole 


i¢€lp5e. 3t.127).39f. (Gal. atereD 


35 to 120 


dsunztlin(roanrcol)ecoetf lis; ) V+4S5um 


129 canrinue 
1190 canrinue 


try (n,n) =3Sun 


100 Contirnue 


90 ¢anrirue 


CALL SJBWOITINE TO [VERSE FILTES FILTERE OSI MSCeRaRi es 


call nonceltrvein) 


86 


lows 


iad 


Nnnonnnn 


NONCAUSAL (CONT. ) 


ADD TH= “ESN OF THE JTNPUT IMAGE TN THE FILTERED RESULT 
45 1435 itO,esize 
do 155 jeO0ecsize 
IrnCisjpdstrylie sp ienean 
155 contiryue 
145 continue 
CREATE AN [MAGE DATA FILE 
cal! suzintfile(in) 
SCALE THE WESULTINS IMAGE ARRAY ANN CONVERT TO BYTE FORM 
call scalelin,integemax,min) 
call intbytelCintessr.in) 
ARTTE THE SENERATED IMAGE INTO A FILE 
O05 150 j=0-rsize 
35 169 3=9,e¢S5ize 
alps deoin( isi) 
1690 contiaue 
wreite(c'i+!) Caln)den=Iersize) 
150 ¢conrirue 
CEdo= FILES 
Close(usaitst) 
close(urit=e2) 
close(3) 
Sr55 
end 


CONV 


THES P2OSR44 GENERATES AV PMAGE TEXTyRE YSING 4 SOUCAUSAL FI2 FILTER 
eee Ce Peston COVALYED F294 A JATe =F [LES TRE FILTER IS 
Doe eA iM biok. SJS2OJTINES SS4LE, ThISYTE, AND SYTEINT ARE 
JSTD TD PREPARE PwWAGE PATA ARRAYS 592 NISPLAY, 


pee TYE VAR TASLES 
Sytceeato: ict )-ore(0tle?,0<)e7) 
rateger nersize,csize,i,jerow,cal,rinters(G:127,dslery),tsrzerntsize 
real*S s~sun,rinl(Or127,02:1le7)-coeflellill,-ll:ll)d,ansunerean,3u79(02] 
ered le7) 
Seer [CES 
overn(unitsti,nares*trarnnaran ral arinsdc. tat isrvoes'aglt',access='dire 
Scr, recoarts12e=352,ragrrer=}{2%) 
overal(unitsed,name='(ratnnaran.irlarinso sat ' rt voe='new’,access='Dire 
‘er’ srecorasize=3e,raxrrer=1e%) 
soen(s-files'trathnrans.gatalconvel.daat',status='ola') 
NESE PARaMETERS 
esize=l27 
csi ze=127 
fsize=1} 
rtSize=elefsize 
READ FILTE2 PaQavETERS J'T) an agray 
On 10 i2fsizer,fsize : 
35 20 json*sizer,fsize 
rea3g(3,25) soef(i,j) 
25 for net (Cic .fe) 
29 cantisaue 
10 caontirue 
S9EAO [¥EGSE TD RE FILTEFEDS INT) At B2RIAyY aN) COVERT TC INTEGER 
45 30 129,ers5i17e 
Peaati; il) Catal), n=0, 127) 
Sa, doy 7 =°.25) 7 
bBinly,j; eal ;,) 
od Conti sue 
390 Contiaue 
Call osyreitanriginrineteeg) 
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CONV (CONT.) 


¢ CONVERT THE INTESE a22 ay TNID & REAL ARQay AND COMPUTE THE WEAN 


™surr}.909 
35 59 120, ersize 
M5 90 4s0,¢s51 76 
rrlie jp Ffloatfintea(1,;)) 
TSUN=In(C 1, ) ONSUM 
50 cantiaue 
S50 canrirue 
Tearsnsunm/float((rsrzetilel(csizetl)) 
erite(*,S3)nean 
SS fornat(420.12) 
c SYATIWATT THE MFA FRM THE ['48GE AaRRay 
35 70 120,esize 
G5 99 j;2U0,¢esrze 
rmnC po zplsrnly, jp demear 
RY CANT AUe 
70 contrrue 
¢ 6F TLCMTER UTE TP AAGE AQRAY 
35 70 nezO,ersize 
G5 103 »=),2size 
Asuel70.2900 
SD 11M tionfsize,fsize 
O> 129 j=nfsize,fsize 
PF 3oasnNe-) 
colame; 
if € (Crow let .0).0r.(ea1.1¢%.9)) aa ts 120 
1¢(Ceow.at.i¢c7).0¢ esl ot le?) son cone 
Asuns(inlronws,colidecaet(i1,j)) ¢45uU" 


120 continue 

110 conrinue 
Sural(aemndtdsun 

100 contiaue 


20 canrtraue 
c SCAL= [at FESULTING [ 44GF ARQ4y AND CONVERT FD SYi aero 
call sca3le(tuno,integ) 
call snrovtelintearain) 
¢ anrmifl= FH= GEVYEVSTED [42G= INITIO A Fite 
d5 157 4=0-rsize , 
35 15) 1=)e 2s ze 
aly )=sin(r,)) 
150 cont sue 
arite(etiel) (Caln)erzdersize) 
150 ¢snti sage 
Cu C52 “lees 
closel(uatr=t) 
c#toselusit=e2) 
closets) 
stood 
ead 


88 


ano an 


SPECOR2 


THIS PROGRAM SCLYES EQUATIONS FOR SPECTRAL CONTENT AND CORRELATION 
IN ONE DIRECTION OF A GIVEN AUTOREGRESSIVE IMAGE MODEL. IT WRITES 
THESE RESULTS TO DEVICES 3 AND < RESPECTIVELY. 

DEFINE VARIABLES 

INTEGER I.K 
REAL“8 PI, ALPHA, THETA. W,ALPHAS .SXW.RXK(-69:469).2-ALPHAK 
DEFINE PARAMETERS 
P1aS.1461592654 
ALPHA#O 999999 
ALPHAS sALPHaaa2 
START X AXIS LOOP AND DEFINE X VALUES 
DO 10 190.99 
We(-1,.08PT)*(2.0MPITa(FLOATC(II/93.0)) 
DO SPECTRAL ANALYSIS 
SXWal 0/01 .0-(2.08ALPHASCOS (W) )*ALPHAS ) 
WRITE(S.15)W.SXw 

18 FORMAT(F10.5.1X.-F10.$) 

10 CONTINUE 

DC CORRELATION FUNCTION SOLUTION 

DO 290 K#0.49 
ALPHAK # ALPHA aa 
RXK(K ) = (ALPHAK/ (1 .0-ALPHAS ) ) 
Ta~1*K 
RXK(T)#RXK(K) 
29 CONTINUE 
CO 30 Ka-69.69 
ZeFLOAT(K) 
WRITE(S.25 2 -RXK(K) 
25 FORMAT(F10.0-.1X.F50.15) 
30 CONTINUE 
STCP 
END 
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ona non 


eno 07 O00 HR HAHAHA 


TeECoks 


THIS PROGRAM SOLVES EQUATIONS FOR SPECTRAL CONTENT AND CORRELATION 
IN ONE DIRECTION OF & GIVEN AUTOREGRESSIVE IMAGE MODEL. IT WRITES 
THESE RESULTS TO DEVICES 3 AND « RESPECTIVELY. 

DEFINE VARIABLES 


INTEGER I1.J-M 
REAL®S P].ALPHA,THETA-ALPHAS.ALSPA.ALPHAS. ALPHAG. THOTN.COS2T We ALB 


®,COS2W. SXW1 »SXW.K eALPHAK. KTHETA,ALPHIK »KP2TMP -KTMP. INTERM.RXK(-49: 
"49).2 
DEFINE PARAMETERS 
PTa3$.14615926586 
ALPHA=0O. $ 
THETAS0.08(PI/12.0) 
IF (THETA.EQ.0.0) THETA#0.000001 
ALPHAS #ALPHAsea§ 
AL SPASALPHAS*ALPHA 
ALPHAS sALPHA#® a2 
ALPHAG sALPHA# 8G 
COS2TeCOS(2.0"8 THETA} 
TWOTHES .O" THETA 
START X AXIS LOOP AND DEFINE XK VALUES 
DO 10 120,99 
Wa(=1.09P])°(2.08P] @(FLOAT(I )/99.0)) 
DO SPECTRAL ANALYSIS 
Aa THET A~wW 
BaTHETAcwW 
COS2w2COS(2.08W) 
SXWl SALPHAS— (COS(A)#ALSPA)* (ALPHAS#COS2T )~(COS(B)*#®AL SPA) *(ALPHAS 
#aCCs2n) 
Swe). 0/01 .O°ALPHAG*(2.08SXwW1 )) 
WRITECS-1S IW SXW 
1s FORMAT(F10.$.1X.F 10.8) 
10 CONTINUE ~ 
OO CORRELATION FUNCTION SOLUTION 
DO $0 420.49 
KsFLOAT(J) 
ALPHAK SALPHASa J 
KTHETAsK® THETA 
KP2TMPa( (2.0°K) 8 THETA )—=PI 
KTMP sK THE TAP] 
INTERMa( COS (KPOTMP )= (ALPHAS#COS (KTMP)))/01.0°ALPHAC-(2.08 ALPHA 
=S*COS(TWOTM))) 
RKK (J) a (ALPHAK/(2.08(SIN( THETA }#82)))"9( (COS(KTHETA)/(1.0—-ALPHA 
®S))¢INTERM) 


GO TO $6 
2S RXK( 5) =e (CALPHAK®8(2.0°K))/02.08(01.0- ALP HAS ) 882) ) 
36 Me-}aJ 


RXK (MP eRKK( J) 
$0 CONTINUE 
DO 40 ya-49,69 
acfuSAT ts? 
WRITES +6S)ORXK(J) 
SS FORMAT(FIC.1+1X.F 30.15) 
40 CONTINUES 
STOP 
END 
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SPECOR3A 


C TH1S PROGRAM SOLVES EQUATIONS FOR SPECTRAL CONTENT AND CORRELATION 
c IN ONE DIRECTION OF A GIVEN AUTOREGRESSIVE IMAGE MODEL. IT WRITES 
c THESE RESULTS TO DEVICES 3 AND & RESPECTIVELY. 
C DEFINE VARIABLES 
INTEGER I,J3-M 
REAL"8 PI.A.8-AB-AMBS.AS.BS, AK] »AK2.BK1] -BKOs INTERM-RXK(~69:69),Z.A 
aK ,KP2.KP] OPAS.OMAS 
C DEFINE PARAMETERS 
PIaS$.141892654 
Az-0.9 
Be-0 .89999 
ABsAuB 
AMBS #(A-B)# #2 
ASsAun? 
BSeBaue2 
C START X AXIS LOOP AND DEFINE X VALUES 
DO 10 180.99 
ba (~1],00PT *(2.0"P IT aC FLOAT(I)/99.0)) 
C 90 SPECTRAL ANALYSIS 
COSlIWeCOS(2.08nW) 
SXW121].0-(2.08(AB*(ASSB)*(CASBS) DECOS(W))°(2,.08ANBSCOSIW) 
SXWe] .0/ 0S KA] ASO (2, 08ANB)S3S°(AS#BS)) 
WRRITECS 1S IW. SXW 
18 FORMAT(F10.$-I1MK-FI0.5) 
19 CONTINUE 
C 3d CORRELATION FUNCTION SOLUTION 
oo 59 Je0.49 
TF(A.ED.8) GO TO S$ 
AK] saue( Jeo) ) 
AKo eAua( yor) 
BK] sBas(J*]) 
BKIsBet( Jeo) ; 
INTERM8(AKI/(1.0-AS )-(( (AK1 28) °ASBK1 2/01 .0~AB) )°CBKI/(01.0-B3) 
*)} 
RXK(J)8(1.0/AMBS) "© INTERM 
GO TO 36 
35 AXeauay) 
KP2sFLOAT(J*2) 
KP 1sFLOAT(J¢1) 
OPAS#1.0°AS 
OMAS#] .0-AS 
RXK (J) 8 ( CRP SEKP 1 PAK) / (02.08 OMAS) DC CFLOAT CU) OKP 1 SAK) / 62.98 OPRAS) 
") 
36 Me~-)]#J 
RXK(M)ERXKC J) 
30 CONTINUE 
DO 4S 3-49.49 
27FLOAT(CJ) 
WRITE(4-65)2.RXK( 5) 
4$ FORMAT(F10.1-1X-F30.15) 
40 CONTINUE 
sToe 
END 


oi. 


YONA nn or tt 


AnnaaAannan 


PIECE 


T47S PIDSAAM TAKER A LARCTOA PIESE O11 ME 2 5} 20S 1 4aGe wooo n ee 
AND PLATES THAT DATS [70a Ae seGe eee 


NETINE varRPAALES 
Ovte 2792591 1),0rm(92511 92:51) ),a0(92127) 
OPEN FE Ees 
osenlunit=l,nagme="lratnnaranvindtrees.tat'-tvoe='old’,access=' aire 
scr’, recordasize=32,ngaunrer=1248) 
Overnlunitsedsnames'(leratnnaranwianltreesbu.tat’»-t(vo9Ve='new's,accesS= dir 
Pest es rerorasire=ztsenanrecrou) 
READ IWMASE TD SE FILTERED [TO 8N BRIAY 
35 30 120,127 
reat(t 'ie*t) Caln),n20,127) 
g5 99 420,127 : 
binfi,jdealj) 
10 cantinue 
39 cantrru.e 
WRITE THE DESIRED Yas 
qd5 tS. 1=55,12e9 
G5 16) j4263,126 
k 5p 7"63 
aa(cdFdIinl(i,;) 


SO2T;O% INTO & FIL=e 


1) 


150 cortiaue 


werte(2'r12*02) (aaln)-21=3,63) 


1S0° € sev 70 


CLOQS ee 2ee5 
closefurnitr=l) 
closetustts?) 
$to2 
end 


INTFILE 


THTS PROSGRav CPESTES & 1284128 OR B4XO4 TUITEGES Plie eae esx pen 
INPUT IMaGE FILE 


DESIYE VER PASCES 
Oyts SCI: 127), 01m 00. 127, a. 127) 
tanteger VepeAs weTerntealls127,02127),sursriian, low 
real*? nean 
OP sr les 
osen(unitzli,nane=' fearhnarn,int tifar2.dat'stvoez'old’,access='d 
sirect srecordasize=se,narrec=t28) 
aoen(2,fiie=’ leatanarn.atatal tinwat dat sstatuss'new') 
*fA0 [¥S5E 7) Se FILISRED INTO Ce aR aay 
35 30 120,127 
reat¢ 116) ata tn).om=0, 127) 
a5. 10°) 720, te7 
o1rnli,;)2a€3) 
40 Cantisaue 
39 cantirue 
CONVERT THE @YTE AQRAY INT) AN INTEGSQ2 amQay 
call avreinr (ring irate) 
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aAanMNrMA NAHM MA ANAANnA AMAA HAHAHA AnAAnHAAAnHA ANA A 


INTFILE (CONT.) 


COMP ITE 415 AND LOW VALJES OF THE 14A4GE DATA AND ARTTE TO TERMINAL 
hiaqas) 
lowe255 
Ae ee Oy ker 
a5 32 j29,i27 
VeCinrealisjireltelowe) lowzrnteals,j) 
PFCingcealeys j)ealeXIIN) NIAN=ZIrAtealr,)) 
32 cContiaue 
31 contre 
wevce(e, 33). HIgh, low 
33 forwatlrr9) 
COvP JTS YVEAN DF INTEGER ARQZAY AND SURTIACT IT FRO THE Data 
suns) 
10 35 1270-127 
Jo 435 329,127 
sur=suntinrealis,;) 
a5 edsr71t nue 
35 soatinge 
nearz=fisar(sund7(129.98129.0) 
aertel(e, $4) 5 yr, 72 aN 
34u fornate(110,e10.3) 
39 30 i2°5,-127 
40 45 F=9,127 
mcerlisfdtintealiejdorvne (reanel.5) 
Us cortinue 
35 sont inye 
mie TNIE SEs €Re4y INTO & 1294128 Cats FILE 
k=12% 
t=123 
arite(2,39)k,! 
32 forngat(2c2i5) 
HO See ies er 
45°60 jevsll2sts 
atytel(e2e TD) 15tealr,pdernrcerlicith),-rateali,p,pe2)eremereai(ry,j+3de3 
enceair,se4)erecealisi ta sinteaeg(i,peddornreanlirejs tT), rtes(r1.,+¢3)019 
mreaglre ft ,irntealier ed M,tntagl(r,ytbido,ryatealry,jel2),rnteaag(s,ye13) 
eo at eg( i, 3! 2) e9rnC eal isi e153) 
on tornat(io15) 
69 Sort rnue 
39 saat Inge 
Gos. -1.f3 
closelusnrt=1) 
close(uratt=2) 
St9o5 
ead 


.e 


TRANS 


TRIS PIOCRAM READS MATS BUT BE A FIL IN FRES FORAST AND CONVERTS 
IT 13 & OPRDNGRAM READABLE FQR 44T, 
integer i,j) 
real? 3(15) 
open(lefrles*fratnnana.gatalagwati. tat ',status='olg') 
oner(2retile=*(ratrnanrnwtaraltfeiwatetar',Starus='new' ) 
Cearsl(let)isjoalld-al2deal3).,al4),a(S)r,alo).al(7),9a(8),a(9),a(10),a( 
@tid,al€l2),a(13),a(14),a(135),a(16) 
d5 10 x«=1,16 
werite(2,20)a(k) 
2° faenar(a29.12) 
10 continue 
$to2 
ead 
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Ana nAnnan 


Q 


NSHP 


THTS P2USREM GENERATES Ay T44GE TEXT IRE "ISIN > Va, TEON ot Se aS eae 
TO AN SUTORSEGRESSISE FILTE? wiISE 2Pa2aDAT TERS AQE HHTAICIED F204 THe 
FILE FILTOE®. SYUSRIITIVE OS4NSS 1S USE> TD Satie aire oe WHITE 
NOTSE AND SUSRIUTINES SCALE AND TNIGYIc @P= JSD Ieee ea eee 
947A ARRAYS FOR OISPLAY. 


JESTINE vaRyJAALes 
ovts a(Jt1e7) sore (0+ ler. ler 
integer Ne SeererSiZercsivZeeonlewMicivie iOleorOwe crs lsinreagl(0:127,03] 
‘27),nom] 
realt3 vall,val2san(Osl27,02127] de tsum,in(Or127-92127) ,coet(-2:9,-2 
139) 

OPEN FILES 
Oner(uNit=l,name=*lratnrarin.,datal en dat set voe='nenw*,access='direct 
:*,recorawsize=se,marrec=123) 
osenl(unitsed,name='lrathrnaran.datalrotate.dat'-type='new',access='3j 
‘rect’ ,racorisi2e=3e,maxrer=12F) 
onernl(setile=*fratnnann.datalfcsa.sat',.-status='314') 
Orern(4e tiles '(ratnrnann.tsaralts.tat',status='new') 

DEFINE PER eet eas 
sees=1234567 
esize=127 
csizesi27 
orl=e 
an1i=4 
nonrleerleorl 

CREATE wALTE NUESE ARRAY 
O> 10 wt2r0,rsize 

ap 29 j20,csrzen-l 
jolrs el 
call saausslse*te-vallevale) 
erlie;)avall 
e(iej,;al)evale 

i) caxtiaue 

10> e5nt 1 ss 

SC sti5 493487 #7 ANO CAVVesT §) 8yte eae 
Call seagtelan,intea) 
call wnerovetelinteraseain) 

ANTTSe THe weETE SMOES= IMAG=e “ARGS y TilieA ee 
95 30 wtOersize 

35 40 j=G,ecsize 
aly)tForn( ry, 5) 
ud eaqntiaue 
erifeflfie4el) falr).r29,eS120) 

30 canrtirue 

2eAD FILTE? PAPAvETEIS FNID 4% ARRAY 
35 30 venonle an 

0% 50 sO, ae! 
reat( 3,55) coef(1,3) 

SS forenat(azc.i2d) 

hu cCantinue 

39 c99ttraue 

APOLy sels SIS? TO 9a. = JIA G ESS eee oa 
49 112 we5,ces8i17z8 

05 12% nz0ersize 
gsun-9.85N 
39.130 steronl, on! 
95 140 7S, 901 

afi (r1.eae3).and.(j.ec.3°)} 95 ts 140 
roOw=4- 1 
colaney 
1¢(Ceagectt.-9).0F. Ceol cite 0) foots lo 
Sdsunmzl(crsetlrepdernlrowscal iy t-i1s4"7 


«fy 
0 


140 c3otinue 

130 eanrinue 
ra(aende(-1.0d08I9¢unN) eanl nem) 

129 ca7tiaue 


110 eantirue 
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m 7 


NSHP (CONT. ) 


SCALS THE 2ESULTINS IMAGE AR24y ANN CONVERT TO BYTE FORM 
cat) scaxlel(in,rnteq) 
call intoveelinrerarsin) 
WRITE THE GENERATED 14A4GE TIO A FILE 
35 159 i=%-rsize 
G5 160 jyF0,esi7e 
aljyltavnli,;) 
1o9 continue 
werce(etiel) (alndss=0-esi12e) 
150 conttrw29e 
GUS: -1cE3 
close(usit=}) 
cCloselunit=2) 
close(3) 
Close(4) 
stoo 
ead 


PGAUSS 


RIVITINE 9549SS(0<,71,22) 

Lar Dowell) 

CSE aa « 
S927TC-2.900 «@ ALQGCRANCK))) 
5.28319590 * RAV(K) 

Soret 5 (3) 

ae 51 405) 


“38 TU (& 


= 
4 
T 


PY ar ee GN Ww mm = AY 
Q a te 
AJ 
- 
a 


<mMmtumy— 
wom 


SCALE 


SyORDUITIANG® SCaxztel sere INnredr,raee WIN) 
iPeweotsa yht we SCULES AY AVeSy TU [*!TeSE Vaud 
vox b4y4 aN) 4g spose 
QESIMNE vaepaaLtes 
mateger racea(J rl 97, D2 1271616 j 
reale?® gre (Ss lLe27e,I2te7T)enriane lows NEAN, Rare DIN 
eee 4) SC atl es 
mhyanse19000,940 
V5w2199090.0 10 
go. 10 1 =0,127 ; 
a3 20 j=6,127 
Mttaertivii.!*@cloe) loetarr(1,}) 
if(agrrli,js).aewnrgn) aranzsarer(y,)) 
sunmsare(is,j)¢sun 
20 Contirue 
10 cantirue 
Mearn=sur4(i 24 .0¢1238.,0) 
arirea(e, 25)5195,]la0e,near 
2S tornatl' oressalet image jsarrae APINE', JI2.5,' low=',dtc.5,' mean= 
ett 26D) 
CUUTTMJE STALIUIG and CONVEST (9 InTESE2 FOR 
359 649«11s0,127 
go 20) 725, 127 
iates(r,-pltrarll (Care lisjslelowsye(CaaxvenialsZ(niqnetow)))+0.500) 
Pytsolt, | y=i1ncea(isi dere (ain) 
40 cContiraue 
30 cantirue 


qr 
iy) 
L 
fn 
—4 

@< 
fry 
fn 
* 

f 

~ 
_ 

< 
rv 
o! 


* 


25) 


Oa A 


SCALE (CONT.) 


CALC ILATE 4164, LOw, AtID MEAW CF SCALED {Mace 
$5,9729.0 40 
Me ansel 1000 .079 
b5e=1I070,040 
a5 30 120,127 
35> 50 7s0,127 
if (loan (inrea(y. |) ional HvaAR TF loarlinreagli,j)) 
rfCfloar(intea(y,3)).1¢. load low=floarclinreal(isj)) 
S4UmM=float(intera(is,s)) esyum 
60 CoONttsuUe 
50 contirnue 
mear=sun/(123.9¢1599,0) 
weice(*, 535) nian, lo5a,meay 
65 torenat(’ sealed inane Jatyz- htgns’. 412.56 Lows", 412.5¢° neem] eel 
eo) 
Petourfrnr 
ent 


INTBYTE 


SYDPFOULINe IntodOvtel(inee#a,a1~) 
T4AYS SISZOFTINE TAsES Aa% [TEESE aORAy A2ND CINVEPTS I[T Ist A 
AXYTE ASRay ; 
DETIVE VAR TESLES 
mrtegqer sejporatea(sas1le7,si27)en 
Sve os nl ele? wore?) 
DEQFIQV FOVVERSI SN 
n=0 
45 10 eo, be? 
95 20 326,127 
TFC Crinteg(rs,jd ett. Gear. Cineegliss).9t.255)) nent] 
iffentes(is 1). ltet) oa 201) =-174 
it (rn tea (1 je ait 25) On Gl es 
if((Cinteali,j).le.t?7) .ang.f(15089(1,);).22.0)) (ones, = 


i ategCi fh) ; eee 
i#(Cineeali,i).9¢.127) <endy enue tp enleee ern bim(i,j)= 


Sinteal1, 4) 2255 

20 cantisue 

1G conetrue 
weice(*, 359) vn 

30 farnar(’ THE NUYSER OF FOTTITS Ju? OF RANGE 15$'°,45) 
retuer 
ead 


96 


aon sa 


nm 


aA nn 


a] 


ANANRANNANnRANNHNA AA 


BYTEINT 


Sydroautine anvreint (arr, inte) 
THIS P2R05R 4M TAKES A AvT= BRAY AND SOSVEST Ses il lowe. THTEGE? ARRAY 


9€71INE VARJASLES 
tatege\er tneea(AsLe7 Ose derejen 
Ovte srn(925127¢931e7) 
PEQFIRY TOVVERSTION 
n=0 
j5. 10 1e0, 2 e7 
dg 29 j=0,12?7 
if€fowmr(1,3)-1¢%.2°123).0re(dim(i,j).atel27)) nsnel 
if€zsin€lr,j)lte.-12"8) jintrealins~jJe-1e8 
1t€sin(iei).at.ie7) isrealiss)=ale7 
it((oinl(is,sleqe.tI1C&®) candi (nin(ises) .1t.0)) wnrealisjde 
‘oinlisj)) #250 
mretown(), i )< 90.0) cand. (Oivel, sielecte?l)) snteali, js )eo1 nls, 5) 
20 contirnue 
10 contirvue 
write(#, 30)n 
30 ¢srnat(' THE YWU4AE2 NF SOTNTS SyT OF RANGE I[5°,15) 
ety s 
eat 


SUBINTFILE 


Suydrostine suarrttriiefis) 
Poses eo il deere AI Dew fete ee? fufse Gta oats He [Le FOQM FILTEIe 4 


JMsce ATA 


Pee tye evan PASE 5 
Pea heen tnt eo Os bed Ne 127), sum 
Pedal *4 wean, i acl le7,y li lel). vane lo6 
eoeN F [LES 
o0exlisfile="fratnnarn. tataldiannar.tatr'.,status=' 
eee gt= 4}ed AH) Lon VSL ylS JF Tae 14aG& Dats aVM a713= TO TE2¥E Nat 
hian=-1,645 
Paws] .94% 
Saees) Neue le?t 
G5 se jecg ler 
TEGCUNC Nile lee ton) loss1r 1, )) 
iC NON stds Qt. 1 an) Aigner nt th, 7) 
32 cantisaue 
31 cantisaue 
werrele, 33) nr an, lode 


55 faenat (2u°6.12) 
Suet 4EGN OF fNTESE2 S#7Avy AND SUZTI4CT ITT F8e4 THE vata 


Now!) 


suns) 
S/S, 35 elle lie 7 
ASS p= dy let 
SUR=suUatINneegli, |; ) 
3) sortinue 
5 roaatinge 
meartflaar(surd/(124.0¢129.9) 
aritele, 34) 5 ym, mean 
343 tornat(ild,e1ad.3) 
JO 35 1=9 Ie? 
30 45 jS,127 
mrealrie,jsleinteali,j yrrne(neane0.5) 
45 coctinue 
35 sortinive 


4 
> 


o7 


AAA N 


” 


SUBTLE Tce 


CONVERT imMaGE nara [YT INTESE2 FORM 
S559 7120127 
Gs 55.) 20 e7 
trtealr, fj dzrne(anlsr,j 29.50) 
$5 CONC iAUe 
3S conrirue 
ARTIS [TNIESER ARQAY TNTO 2 {248x128 OaTa FILE 
xes128 
bs1 28 
wrice(1,39)&,1 
39 tyenae (215) 
35 350 120,127 
33 56 p=0 112,14 
wriredl, T° )inteali,j),integ(ini+tl,imtealss | oe) nee eee 
emcee strc ys tte tatea(t. 145) einteals. jo), int eat 1, je) lems 9 
eCeg(te Fei atea(s1s F410), i ace slid tl), anes alee cee #13) 
eernrermiri tla), rinceglisi #15) ~ 
70 fornat(!o1S) 
60 cot irnue 
50 ¢antirue 
SEIS: eee s 
Close(luntt=]) | 
Pefyrr 
eat 


NONC 


SuDrFOutINe nornc(tCery,1n) 
THIS SJSLO TINE TAKES Aan [4AGE ASREY AND FE TETERS SI Gea ieee 


NONGSUSh FIP a) yen 


DEF VE vaARPAaLles 
Vateger nersrzer,cstze-isg,jrefFow, cnt, fsize,mfsizersiandecl,indere 
Fealt? yssun,rn(I2-127 652127) ,coV#tleloslus 19s lod ensunernean,dumaleols 
9128, =) ites) .niagns lowe try (0<le77O03te7) 
OPE Ares : 
osen(letile=* (ratnnana.sgacalacfiles2l .dat'estatus= olf’) 
DES INE PARAMETERS 
Ader 1S9 
irtenrZ2=j]27 
esizerzi2/ 
Csizer127 
fsi ze 219 
sfsr1zer=eletfsize 
RESD FILTE? SARAMETE®SS INT) Aw ARRAY 
35 10 yvemtsize,fsize 
G5 20 jsmtsizerfsize 
reay(1-25) coef(t,j) 
25 forrat(329.12) 
20 cantisaue 
10 contirue 
TAKE Tee [vO yt AgRAr 48ND CJ42UTC THE AE AN, HIGA, 
asur=).000 
hign=-19000.0730 
low210090.040 
O05 590 121Nndeckls,intere 
Q5 50 jytindenrls intend 
Sunster(ie jy) ¢nSUD 
ifCecevliszp egtehighl Nvanstery(i,)) 
1f{Ceryliei).tt~loa) lowetrv(i,;) 
60 Co Ntirave 


50 caontirxue 
mearsnsym/tloat (Cintenrcdmiateni¢lde(intercd-irndgertel)) 


wricet(*#,55) ean, nrvansloa 
55 fornat(' snout inage arrav= means 
a2.5) 


AND LOW VaLUES 


*,d12.5,° hhoRe [ole v5. lows*,dl 


98 


Cc 


NOnnnann 


NONC (CONT. ) 


FILTER THE [WAGE ARRAY 
35 90 nNz=Indexnlringened 
G9 1990 wneintzexnl,tndex2 
Osun=0.070 
G> 110 s:eHfsiz7e,fsize 
d5 120 j=anfsize,fsize 
Powst Ne) 
Jone, 
Pie ee erect a) Conagextl S500 1e0 
VEC (rows gt eirexc«e).orw(cesteat.intex«2)) 30 to 120 
Osur=(tervlrow, cal) ecoeflr,j)) +36un 


120 continue 
110 continue 

ir(Casendrotsun 
100 Contisaue 


‘ e 
Memes oa: 139X130 aRRAY INTO A 128x128 AVRAY, JF NECESSARY 
jo 130 i20,127 
Jo 140 j=0,127 
wa(i,plesgunol(i,j) 
149 corztinuve 
130 continge 
Close(}) 
retut 
end 
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MAKMAT 


SRR OY OG aes tc ena 
J borne 
Aj “3 41 66 60 


ee tA li, aes ee 


ALO: ees 
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APPENDIX B 
DERIVATION OF THE POWER SPECTRUM AND AUTOCORRELATION 
FUNCTION FOR THE TWO POLE AUTOREGRESSIVE MODEL 


Power Spectrum 


H(eJ”)H(e7J%) 


IH(eJ®) |2 = Hed”) -(H(eJ”))* 


In this case: 


1 ‘le 
Hi(Go a c=, H(e73 9) = ee 
1 +ae be 1 +ae3” 
Calculating H(e!” )H(e7J%): 
oe Eee 
(itae 2%) (l4ae2”) Teel gee ee pe coane ee 
The final result is: 
: a 
Sy (5) = SIH(e29) |= = See (B.1) 
(l+a~)+2acos (w) 
Autocorrelation Function 
Starting with H(z) for this case: 
1 +az 
Per Ref. 7:p. 158; 
Z-1(H(z)] = h(n) = (-a)P-u(n) for a<1 


(u(n) is the unit step function) 


Per Ref. 5:pp. 391-395, for the white noise input case: 


Ry(2) = > oh i) sini (B.2) 


n=-© 


Substituting h(n) above into Eq. B.2 
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CoO 





n n= 
R(t) = ) (a) = (-a) & >? 0 
i) ee 
(-a)~ n= 


The summation term may also be expressed as: 


co 0 va 
y (-c1) 2n = y (20) 2n = i (-c1) en 


Per Ref.12:p. 8, the summation terms on the right are equal 





Os 
co Q-1 22 
2n i! 2n 1 - (-a) i 
(-a) = > «Sa <1) Y (-0) = (a < 1) 
ade 1 = (-a) n=0 1 - (-a) 
As a result: 
00 2x eae 
y (=a) 29 = —— gee?) Pe aL as ea) 
ae eC) ae) 1 
Substituting and using (-a)* = a2 yields: 
2 ae 
Ry(2) = - ae = Je ee Oe anced y= il (B.3) 
(-a) 1 -a 1 -a 


NOES 


APPENDIX C 


GRAPHICAL RESULTS FOR THE POWER SPECTRUM AND AUTOCORRELATION 
FUNCTION FOR THE TWO POLE AUTOREGRESSIVE MODEL 


Power Spectrum 





ALPHA=0.95 ALPHA=O%3 
© 
(2) rs) 
2 o 
Eg Ee 
Zz Fed 
ts 3 
< =< 
= > 
o o 
=3.16 0.00 3.15 -3.16 0.00 3.15 
OMEGA OMEGA 
ALPHA=—0.95 ALPHA=Os8 
oS 
ee tx) 
= = 
= = 2 
a = 
= = 
= x 
© 
= 15 3.16 -3.16 0.00 3.55 
OMEGA | OMEGA 
ALPHA=0:25 ALPHA=—0.25 
: a 
ral a 
= =5 
: E- | 
Bea oe 
O C5 
< < 
2 < | 
© ; 
-3.18 0.00 -3.15 0.06 3.15 
OMEGA OMEGA 


GG 


ALPHA=0.75 ALPHA=0.1 


30 
z 


MAGNITUDE 
16 

MAGNITUDE 
{ 


0 
Q 


-3.16 0.00 3.15 -3.16 0.00 3.16 
OMEGA OMEGA 


Autocorrelation Function 


ALPHA=0.5 ALPHA=0.95 





" 60 0 60 


L 
at A=—O.9D 





107 


ALPHA=0.25 ALPHA=-0.25 
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APPENDIX D 


DERIVATION OF THE POWER SPECTRUM AND AUTOCORRELATION 
FUNCTION FOR THE FOUR POLE AUTOREGRESSIVE MODEL 





Power Spectrum 


IH(eJ) |2 = H(eJ%) -H(e7v) 
Tm this case: 
; il 
ae) = -76 = ojo.  . ..3(Ot»)  —-jJ(O-w). 2 j2w 
(l-ae J ae miGi=aee en) l-ae —1e +a e 
: 1 
ages) = =o = loa = —j (8+w) 7 (8-—w) , 2 -j2w 
(l-ae J¥%e J”) (1-get%e I) l-oe -ce +0,e 


Multiplying the above expressions yields: 


Ih 


L-aed PP) Loe JOH) Ae eee eue a eae ae 


eA) 2, 2 ~J20_. 3-9 (0-w) 


= ged (9-4) Oe 


aL 


n tes eg red (8+W)_ 3,75 (8+) 4 


Seombining terms: 


Hl 


Me)%)-H(e 2) = ) , ee 
J (6-w) 43 (8-W) ) 442 (92947528) 


eee ie 


_ on (ed Ete) 7 (OW) 5 12 (JAH, 20 


~ (ed (OH) 47 (8-W) ) _ 9 (3 (OH) 75 (8+u) | 


Using Euler's relation and combining terms: 


al 


Ne a a Re A 
1l+q “+2 (a -cos (6-w) [a~ +a] +a~cos (268) -cos (6+w) 


x ee eaeos (oun 


OD 


Using cos(8-w) + cos(8+w) = 2cos(8)cos(w), and Sines 
Sy () = H(e)”) -H(e7I%) + 02, with o@ = 1 the final result is: 


Sy(w) = eee 


1t0*+2 (a*—2 [a--0}cos (8) cos (w) 40~ (cos (26) coe ane 


Autocorrelation Function 


yee eee a eee 1-(e2 84609) 2 22? 
Using Euler's relation: 


1 
H(Z)) = ———>—_>—,— [Ze are 
1-2acos (8) z Tez 2 


Per Ref. 7:pp. 204-216, partial fraction expansion can be 

used to find the inverse Z transform. To do so H(z) can be 
expressed in the form: 

ar Zz 

H(Z) = ——_————~— = 

£2) @ -l 38 


(oqems oa diesel oem (enema eae. 


Using the partial fraction expansion and table look up [Ref. 


7:p. 158] yields: 


n 
oop eee alle. 
h(n) Sin (6) cos (n6+6 5) U(n)) tore nl (D.2) 
Since cos(6- =) = sin(o9), the final expression for h(n) 1s: 
yn 
= 1 ° Das 
h(n) Sint) sin ( (n+l) 6) -u(n) ome ( ) 


iO 


Boresimplicity in further derivation of Ry (2) based on h(n), 
Eq. D.2 will be used. 

Using the expression for the autocorrelation function of 
a random process represented by the above filter with a 


white noise input [Ref. 4:pp. 391-395]: 


co —2 co 
Ry(2) = J hin)-hineey) = —%— J aM" -cos (note) 
iy ae Z 
n= sin (6) n= 
° cos ( (n-2) 6+6—5) Co.) 
n= % in the summation index since h(n) is causal. 2& is 
Peoumea to be greater than zero here. For %&< 0, Ry(%) = 


Ry (-% ) by symmetry of the autocorrelation function ([Ref. 
5:p. 388], so we can proceed assuming only positive values 
Or 2. 


Using the trigonometric identity for a product of 





cosines: 
ee l 
Ry(2) = —=s— ) of [% -cos(26)+5cos (2n6-26+26-7)] 2 > 0 
uf - 2 2 = 
Sin (8) n=2 
cos (26) «a on of” y 2n 5 
R_(2) SS ee ee aM". (2n6-26+20-7) = (D-9) 
Y 2sin® (9) n= 2sin*(6) n=2 
Lee, 
co co tek | 
Using ) a8 ) = ) and standard geometric progression 
n=2 n=0 n=0 
identities [Ref. ]2:p. 8]: 
Co aa Q 
|- 
) ge — and . a — 
n=0 o n=O 


ie 


For the first term in Ry (2) 











Son On ue: “ae iH ieee 
) a = ) Yoeet=) ) ae 7 Di 2 
n=2 =() n=O l-a 1l-a 
yee | 
wi D.6) 
ana 


For the second term in Ry (2): 
Lek d = —- £6+ 26 -T 


Using Euler's relation: 


j (2n6+¢) —F (2m) 
cos(2n0+)) = = 5 = 


y 07 cos (2n6+4) = a y eae ae i. 73 (2n8+9) | 
n=2% ne 
= = [geen ee 922 g7d (2né+6) (D7) 
n=k n=< 


For large n, it is evident that the a2" term will tend 
to make the term in each sum approach 0 for a< 1, and thus 
ensures convergence and a closed form expression for each 
term. Pursuing the mathematics required to find this closed 
form expression we have: 


val 
a y (eI?) 29) 
n=0 


2 yO . 
: ; ye gd (2né+o) y ee ee 
n= 


eae 
—-[ ) (ae 
g 2 “'n=0 


n= 


ia 


y 2n  3(2n6+o) _ et ib Lote! 
L OF satis or eee, ee 0. 2 
n= 1-(ae-~) 1-(ae 
Pearce el (D.8) 
ee 


For the conjugate term we must have: 


j oan ed (2no+$) eee 


het Dees) 


L 
2 
Next a common denominator must be found to sum these two 


terms: 


92% oJ (28240) ge 
76,2 


yo (1-(ae 


om 
“3°)4) 

ores (282+) _. (2242) 5 (202-2646) 
76, 2_ 


2 (1- (ae 


Me es = (oe yaa) 


oe vee (2Ok+6) Me (aet 


D 38,2, 


yee q7 5 (2024) __ (2242) .—j (282-2049) 
2 


ae tce 2°) ) (1~(ae 


2 Glatere!) \ eine 
By Euler's relation: 
ees cele tas = Tepdeeoe aio: 


Adding the terms with the common denominator yields: 


92% fod (262+) _ 2,5 (202-2844) | -j (202+) _ 2.-j (202-2644), 


Gee caso ac. | 


dele 


Again using Euler's relation the above expression reduces 


tO: 


o*" [cos (262+) -0“cos (202-20+9) 1 (D.10) 


[eonceceree mae 


which is the sum of the last two terms in Eq. (Digie 
Substituting Eq. (D.10) and Eq. (D.6) into Eq. (D.5) yields: 


Bonen nee 
St 


2sin- (6) l-o 





RS) 
4 


o_o [eos (20245) a" cos (202-204) (41) 


+ ——————————_ 
Dea) )2eeeeeeatom 


Combining and canceling terms and substituting for ¢ and 


noticing that the same result must hold for 2 < 0 we have: 


a! _cost|t |e) ,cos((2+|2])0-n)-a“cos({2/8-m), (p12) 


R (2) = 5 
Y 2sin* (8) l-a 1-20,“cos (26) +a 
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APPENDIX E 


GRAPHICAL RESULTS FOR THE POWER SPECTRUM AND AUTOCORRELATION 
FUNCTION FOR THE FOUR POLE AUTOREGRESSIVE MODEL 


Power Spectrum 
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APPENDIX F 


DERIVATION OF THE POWER SPECTRUM AND AUTOCORRELATION 
FUNCTION FOR THE FOUR POLE AUTOREGRESSIVE MODEL 


(WITH TWO POLES ON THE REAL AXIS) 


Power Spectrum 


Using initial results from Appendix D with the necessary 


modifications (including §= 0) we have: 


H(e7J%) = 5, 

1-a,¢€ ce +0, 04, 
: 8 

nce — a 
ai eee Jw 20 
I-a,e ~ -%,e ~ +a,%¢ 

Therefore: 
i 


H(eJ“) -H(e7J®) = 3 3 =p ee 
oI, Ade J4W_ JW 
ak ae ad too, © ae tor to, 


2 25 Jw 2a 
O-, O4,¢ —O4,€ ui mee a4, ¢ 

320.2. jw 0 2 eee 

+ a O° 01% © oe pce 


Combining terms and using Euler's relation: 


ils 
1 (ai, toy 41204, 41, o1-) * 2008 (Ww) +201, 04,608 (213) 415 





H(eJ)-H(e7J%) = 
= 2a ooo ee 
tb Dab 


20 


Assuming o* = 1 and “since Sy (w) = Hes yon(e J yeas. the 


final result is: 




















alk 
a i ie aa >; >. 2 
¥ 1-2 (or, tay, toy, +1, ch.) COS (w) +201, cy COS (2u) +0420, oy, 
eee (ped) 
+ oO) +0. 04, 
Autocorrelation Function 
i 
H(z) = ; ; — 
cme ee oue 2 
1 ae" 2 ape ~zZ +O, 0,2 
Letting 8 = 0 
Z 
iL Z 
H(2)) = ———-+,-——_5 = 5 
1-(a+o,) 2 Hol O42 z — (OL FO}, ) Z+0._ 0}, 
Iz| > [a 
Expanding in terms of partial fractions we have: 
a Z. nt Z 
eee ee 
"arb eee Dec “b 
This corresponds to the impulse response 
on 
2) a - oh) su (n) Ge< 1 (Ee) 
Soe 0 7a 


Proceeding as in Appendix D: 


te 


n+1 n+l yntl-2 n+l—-2 
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Continuing with the same principles and assumptions as in 


Appendix D, we have: 
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Making the appropriate substitutions in the expression for 


Ry (2), we Nave: 
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Combining terms yields the final expression: 
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APPENDIX G 


GRAPHICAL RESULTS FOR THE POWER SPECTRUM AND AUTO- 
CORRELATION FUNCTION FOR THE FOUR POLE AUTO- 
REGRESSIVE MODEL (WITH TWO POLES ON THE REAL AXIS 
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LAPLACIAN INVERSE FILTER FORMS 
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APPENDIX I 


CONVOLUTION OF LAPALACIAN DIFFERENCE OPERATOR AND 
VARIOUS SIZE FIR INVERSE FILTERS 
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